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1 .  Introduction 


The  aim  of  this  work  was  to  develop  efficient  algorithms  and 
programs  for  the  solutions  of  the  ionospheric  diffraction  and 
statistical  radiotomography  (RT)  problems.  The  diffraction  and 
statistical  RT  problems  are  related  with  reconstruction  of  the 
scattering  ionospheric  irregularities  structure.  The  use  of 
tomographic  methods  is  perhaps  an  unavoidable  stage  in  evolution 
of  nearly  all  diagnostic  systems.  At  a  sufficiently  high  level  of 
development  of  remote  sensing  technology  and  data  processing 
resources,  it  become  possible  to  reconstruct  the  spatial  structure 
of  a  medium  on  the  basis  of  tomography.  Tomographic  approaches 
have  already  transformed  sensing  methods  in  many  fields  and 
provided  fundamentally  new  results.  The  major  achievements  of 
tomography  in  medicine  and  molecular  biology  are  widely  known. 
Tomographic  methods  have  made  it  possible  to  detect  previously 
unknown  phenomena  in  geophysics  (seismotomography  of  the  Earth  and 
acoustic  tomography  of  the  ocean).  At  the  present  stage,  radio 
sounding  technology  makes  it  possible  to  use  satellite  resources 
to  perform  ionospheric  sounding  in  a  wide  range  of  different 
positions  of  the  transmitting  and  receiving  systems  and  to  use 
tomographic  methods.  In  connection  with  this,  work  on  RT  of  the 
ionosphere  has  began  actively  in  recent  years. 

The  purpose  of  this  paper  is  to  describe  shortly  the  RT 
methods  El -8]  and  program  developed  for  reconstruction  ionospheric 
irregularities  structures  together  with  some  results  of  computer 
simulation.  We  are  dealing  with  reconstruction  of  the 
irregularities  of  the  electron  density  and  of  the  effective 
collision  frequency  by  scattered  radio  waves  field.  Mathematicians 
used  very  often  designations  the  inverse  problem  (IP)  or  the 


inverse  scattering  problem  (ISP)  for  such  problems  of  structure 
reconstruction  owing  to  the  scattered,  field.  Here  we  used  the 
terms  IP  or  ISP  and  diffraction  tomography  as  synonims.  The 
ionosphere  has  a  rather  complex  structure;  local  irregulariries  of 
various  scales,  including  turbulent  regions,  are  present  in 
addition  to  a  quasi-stratified  background  with  large  scales. 
Therefore,  problems  of  satellite  RT  of  the  ionosphere  should  be 
divided  into  statistical  RT  and  deterministic  RT  problems.  The 
latter,  in  turn,  are  divided  into  problems  of  diffraction  RT  and 
ray  RT  of  large  structures,  where  diffraction  effects  are  not 
significant. 

There  is  no  special  reason  to  provide  a  strict  classification 
of  any  field  of  science  developing,  but  certain  delimitations  and 
explanations  of  the  terminology  must  be  provided  so  that  we  will 
not  later  be  reproached  for  "inventing"  new  terms  unnecessarily. 
In  the  cases  where  any  projections  or  "cross  sections"  of  an 
inhomogeneous  object  (or  a  certain  transformation  of  the  object, 
such  as  Pourier  transform)  are  known  from  remote  sensing  data  and 
the  problem  is  to  reconstruct  the  structure  of  the  object,  the 
problem  should  be  considered  as  a  tomographic  problem.  At  present, 
the  mathematical  foundation  for  tomography  is  related  to  integral 
geometry,  where  it  is  necessary  to  reconstruct  an  object  using 
data  on  it  in  the  form  of  integrals  with  respect  to 
small-dimension  manifolds.  Therefore,  the  term  tomography  is 
understood  not  in  narrow  initial  sense,  as  layer-by-layer  study  of 
the  structure  of  the  inhomogeneous  objects,  but  ruther  in  the 
wider  sense,  as  recording  of  projections  or  cross-sections  of  an 
object  and  subsequent  reconstruction  of  the  structure  of  the 
object  from  them.  The  projections  of  an  object  are  various 
integrals  with  respect  to  small-dimension  manifolds.  The  range  of 


problems  of  radio  remote  sensing  of  the  ionosphere  using 
satellites  examined  in  this  paper  includes  reconstruction  of  the 
structure  of  ionospheric  irregularities  from  files  of  various 
tomographic-type  data  (projections,  cross  sections),  which 
provides  the  basis  for  using  the  term  RT  of  the  ionosphere. 

Problems  of  satellite  RT  should  be  divided  into  deterministic 
and  statistical  problems.  In  the  case  of  a  deterministic  RT 
problem,  it  is  necessary  to  reconstruct  the  structure  of  certain 
large  irregularities  or  a  group  of  irregularities.  If  a  large 
number  of  irregularities  occupies  a  certain  region  in  space,  it  is 
unadvised  to  reconstruct  the  structure  of  all  realizations 
(generations)  of  irregularities  every  moment  of  time.  Here,  it 
makes  sense  to  pose  the  problem  of  reconstruction  of  the 
characteristics  of  full  ensemble  of  irregularities,  the  structure 
of  the  statistical  characteristics  of  the  irregularities  such  as 
the  correlation  function  of  the  electron  density,  etc. 

We  will  single  out  a  number  of  fundamental  features  of 
problems  of  RT  of  the  ionosphere  or  the  IP  of  reconstruction  of 
the  structure  of  the  inhomogeneous  ionosphere.  The  dimensions  of 
the  transmitting  and  receiving  systems  feasible  in  practice  are 
much  less  than  the  distances  from  them  to  the  irregularities  to  be 
reconstructed,  which  are  hundreds  of  kilometers,  i.e.,  the 
aperture  angles  are  small.  Since  it  is  extremely  difficult  and 
expensive  to  create  transmitting  and  receiving  systems  with  a 
large  number  of  receivers,  the  necessity  for  aperture  synthesis  in 
one  coordinate  becomes  clear.  Synthesis  apertures  can  be  realized, 
for  example,  using  a  moving  transmitter  on  a  satellite.  We  will 
emphase  that  satellite  radio  sounding  makes  it  possible  in 
practice  to  obtain  files  of  various  tomographic  data  and  to 
actually  realize  RT  of  the  inhomogeneous  ionosphere.  In  view  of 


the  small  aperture  angles,  it  is  advisable  to  pose  the  problem  of 
reconstruction  of  the  structure  of  irregularities  with  dimensions 
singficantly  exceeding  the  wavelength.  Therefore,  the  IP  of 
reconstruction  of  the  structure  of  scatterers  which  are  large  in 
compasion  with  the  wavelength  will  be  examined  here,  but  the 
dimensions  of  the  irregularities  may  be  both  greater  than  or  less 
than  Fresnel  zone,  and  diffraction  effects  must  be  considered  in  a 
number  of  cases. 


2.  Physical  and  mathematical  formulation  of  the  radio  tomographic 
sounding  problem. 

Radio  wave  propagation  in  near-earth  space  and  scattering  in 
plasma  irregularities  in  the  Ionosphere  are  discribed  by  a  system 
of  Maxwell  equations  with  the  corresponding  material  equations 
[93.  In  the  general  case,  the  dielectric  permittivity  of  the 
plasma  is  a  tensor  value  due  to  spatial  dispersion,  even  In  the 
absence  of  a  magnetic  field.  However,  the  spatial  dispersion  can 
be  disregarded  for  ionospheric  irregularities  of  both  natural  and 
artificial  origin,  since  the  thermal  velocities  of  the  electrons 
in  the  plasma  are  significantly  less  then  than  the  speed  of  light, 
i.e.,  the  approximation  of  a  "cold”  plasma  is  valid.  likewise,  the 
velocities  of  diffusion,  mixing  and  other  transport  processes,  as 
well  as  the  velocity  of  the  transmitting  and  receiving  systems, 
does  not  exceed  10  km/s  (v/c  <  3'10-4),  therefore  ,  the 
quasistationary  approximation  is  valid;  within  the  framework  of 
it,  the  "slow"  time  dependence  of  the  characteristics  of  the 
medium  and  fields  can  be  considered  as  a  dependence  on  the 
parameter. 


In  place  of  the  equations  for  the  field  £(r,t)  of  a  radio 
wave,  it  is  convenient  to  use  equations  for  the  complex  amplitudes 
t  of  the  corresponding  monochromatic  components,  i.e.,  it  is 
advisable  to  introduce 

£(  r,t)  =  S(r,t)e"iCi3t  +  3*(r,t)eia}t. 

For  brevity,  the  complex  amplitudes  S(r,t)  with  a  "slow"  time  will 
be  called  the  field  hereafter.  Subject  to  the  observations  made 
above,  within  the  framework  of  the  quasistationary  approximation 
and  a  "cold"  plasma,  from  the  Maxwell  equations  we  have  the 
equation  for  the  field 

2 

-*  (i)  A  ->  -*■ 

AE  +  — p-  e  E  -  grad  div  E  =  0,  (1  ) 

cr 

where  s  is  the  complex  dielectric  permittivitj^  tensor,  f  =(±i/2%  is 
the  frequency.  Analysis  of  radio  wave  propagation  in  an 

A 

inhomogeneous  magnetoactive  plasma  with  a  tensor  dependence  of  s 

is  a  difficult  problem.  There  are  no  algorithms  for  calculating 

the  parameters  of  radio  signals  and  the  characteristics  of  the 

scattered  fields  in  the  general  case  of  an  arbitrary  dependense  of 

e(r)  on  the  coordinates.  Therefore,  the  statement  of  the  ISP  of 

reconstruction  of  the  structure  of  the  tensor  characterizing  an 

ionospheric  irregularity  using  data  on  the  measured  field  is 

unrealistic  at  present,  and  hardly  advisable  in  general.  The 

non-diagonal  elements  of  the  tensor  can  be  disregarded  at  high 

frequencies,  since  they  do  not  exceed  the  square  of  the  ratio  of 

the  plasma  frequency  fN  to  the  sounding  frequency,  i.e.,  ~ 

(fN/f )2(fH/f )  (fH  is  the  gyrofrequency)  [93.  For  example,  the 

6  -3 

typical  maximum  concentration  in  the  ionosphere  is  NQ  ~  10  cm  . 
For  it  and  where  f  >  50  Mhs,  (fN/f )2(fH/f )  <  0,001.  Likewise,  the 
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last,  "depolarization"  term  of  equation  (1),  the  order  of  which  is 
determined  by  the  ratio  of  the  emission  wavelength  to  the 

characteristic  scale  of  variation  of  the  concentration,  can  also 

be  disregarded  at  high  frequencies. 

Thus,  in  the  case  of  high  sounding  frequencies  vector 

equation  (1 )  decomposes  into  three  scalar  equations,  and  it  is 
sufficient  to  examine  the  equation  for  one  field  component 

A  E  +  k2s(r,k)  E  =  0,  (2) 

N  ( I* ) 

where  e(r,0))  =  1  -  — p - e - ~ - »  E  =  2'jcf/c  =  Q)/c  is  the 

k2(1  +  iv(r)/u>) 

wave  number,  r  is  the  classical  electron  radius,  v(r)  is 
effective  collision  frequency,  N(r)  is  the  electron  concentration; 
"ion"  contribution  to  the  dielectric  permittivity  can  be 
disregarded  [33.  The  relation  for  s  has  the  same  form  in  both  the 
International  System  of  Units  (SI)  and  CGS  system;  only  the 
expression  for  re  changes. 

Then,  it  will  be  convinient  to  introduce  the  complex  function 
q(r,w)  =  k2(1-e)=  4TCi*eN(r)/(1  +  iv(r)/oj)  and  divide  it  into  two 
parts:  q  ->  q0(z,u))  +  q(r,w),  one  of  which  corresponds  to  the 
regular  stratified  ionosphere  with  the  dependencies  NQ(z)  and 
vQ(z),  and  the  other  of  which  is  related  to  the  three-dimension 
(3D)  irregularities  on  the  background  of  the  stratified  medium 
N(r)  and  v(v).  The  specific  and  concrete  form  of  relations  for 
vQ(z)  and  v(r)  depends  on  the  specifics  of  the  problem  [13.  By 
reconstructing  the  structure  of  the  complex  function  q(r,u),  it  is 
also  possible  to  reconstruct  the  structure  of  these  3D 
disturbances  with  known  information  on  the  regular  large  scale 
ionosphere.  Therefore  we  will  deal  with  the  two  problems:  first, 
the  problem  of  reconstruction  of  the  regular  large  scale 
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ionosphere  q(Q ( z ,  w )  -  ray  RT  problem;  second.,  the  problem  of 
reconstruction  the  2D  or  3D  irregularities  structure  - 
diffractional  and  (or)  statistical  RT  problem.  The  introduction  of 
q  is  also  justified  by  the  fact  that  the  quantity  q  ™  (1  -  a) 
applies  to  the  generalised  susceptibilities  [10]  characterising 
the  response  of  the  medium  to  the  field.  Furthermore,  the  scalar 
Helmholts  eguation  assumes  the  form  of  stationary  Schrodinger 
equation,  where  q  is  the  complex  potential: 

AE  +  k2E  -  q0(s,o))E  -  q(r,u)E  =  5(r  -  rQ)  (3) 

The  delta  function  in  the  right  side  of  (3)  corresponds  to  the 
point  sounding  source  in  the  majority  of  satellite  radio  sounding 
experiments,  since  the  wave  can  be  assumed  to  be  spherical  within 
the  major  lobe  of  the  directional  pattern  of  an  arbitrary  source. 
The  dimensional  constant  Aq  in  front  of  the  6-function  in  equation 
(3)  is  omitted  so  that  it  need  not  be  rewritten  in  all  subsequent 
formulas.  The  constant  is  expressed  by  means  of  the  total  power  of 
a  point  source  [13.  Here,  AQ  corresponds  to  the  amplitude  value  of 
the  field  of  a  point  source,  i.e., 

E0(r)  =  -  Aq  exp(ikr)/(47ur). 

Hereafter,  we  can  always  return  to  the  dimensional  field  E  -»  AQE 
when  necessary. 

The  ISP  for  equation  (3)  can  be  formulated  in  the  following 
manner:  it  is  necessary  to  reconstruct  the  structure  of  the 
irregularities  from  measurements  of  the  field  in  a  certain  limited 
region  of  a  fixed  surface  in  a  limited  range  of  frequencies  and 
positions  of  the  sounding  source.  Knowledge  of  the  complex 
function  q(r,co)  makes  it  possible  to  reconstruct  the  structure  of 
both  N(r)  and  v{r).  Here  q(r,w)  is  a  finite  function,  3ince  the 
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irregularities  of  N(r)  are  confined  in  space.  The  regular 
stratified  ionosphere  is  assumed  to  be  known  here,  i.e.,  the  IP  of 
reconstruction  of  the  regular  ionosphere  is  solved  by  other, 
sufficiently  well-developed  methods,  or  its  influence  can  be 
disregarded,  for  example,  as  in  the  UHP/VHF  waveband. 

The  sheme  of  the  experiments  on  diffraction  and  statistical 
RT  is  presented  in  Fig.1.  The  satellite,  with  a  transmitter  on 
board  moves  at  altitude  z  =zQ;  the  receiving  system  is  located  in 
the  plane  z  =  zR  (on  the  ground).  The  receiving  system  can  be  a 
transverse  array  of  receivers  for  diffraction  RT  or  statistical  RT 
(or  a  collection  of  transverse  arrays).  Receivers  of  array 
arranged  along  (the  y-axis)  a  line  transverse  to  the  satellite 
path  (the  x-axis),  arrays  being  separated  by  distances  of  hundreds 
of  kilometers.  The  scatterer  (a  group  of  irregularities)  is 
located  near  the  altitude  z  =  zg. 


3.  The  solution  of  the  diffraction  radiotomography  problem 
(inverse  scattering  problem). 

Rigorous  results  of  solution  of  IP  (3)  for  the  complex 
potential  do  not  exist.  We  will  emphasize  that  in  the  case  of 
propagation  of  waves  of  various  nature  in  the  actual  physical 
media  the  complex  nature  of  the  potential  q  has  a  fundamental 
character  and  the  corresponding  dispersion  relations  relating  the 
real  and  imaginary  parts  of  the  potential  can  be  derived  [103.  The 
fundamental  results  on  solution  of  the  3D  inverse  problem  for  the 
Schrodinger  equation  with  the  real  potential  q(r)  not  depend  on 
the  frequency  (where  q0(r,w)  =  0)  obtained  in  [11,123  should  be 
noted  here.  As  to  the  complex  nature  of  the  potentials,  there  no 
rigorous  results  for  the  ID  IP,  i.e  the  extension  of  the  familiar 
Gelfand-Levitan-Marchenko  algorithm  to  this  case  is  unknown,  even 
more  so  for  a  arbitrary  dependence  of  the  potential  from  the 
energy  (frequency  w)  [133. 

Approximate  approaches  to  solution  of  the  3D  ISP  are 
therefore  of  undoubted  interest.  The  approximate  approach  to 
solution  of  inverse  RT  problem  is  justified  by  the  properties  of 
the  regular  Ionosphere,  the  typical  scales  of  essential 
disturbances  of  which  significantly  exceed  the  radio  sounding 
wavelength.  Therefore,  the  geometrical  optics  approximation  is 
applicable  to  description  of  wave  propagation  on  the  background  of 
the  regular  Ionosphere.  The  scattering  In  weak  irregularities  is 
described  within  the  framework  of  such  familiar  approximations  as 
the  Born  approximation  and  Rytov  approximation  (BA  and  RA).  More 
specialized  methods  have  been  used  to  calculate  the  scattering  in 
strong-scattering  irregularities. 

There  are  many  papers  connected  with  the  approximate 
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approaches  to  ISP  [14-18].  The  scattering  by  weak  Irregularities 
was  described  within  the  BA  and  RA.  The  iterative  methods  were 
used  for  strong-scattering  irregularities.  But  there  are  some 
specific  features  of  the  ionosphere  problems,  namely:  the  size  of 
irregularities  is  very  big  (it  can  exceed  the  probing  wavelength 
in  IQ3-  104  times),  the  ionosphere  irregularities  can  be  either 
weak-scattering  or  strong-scattering.  It  is  very  difficult  to  use 
iterative  methods  for  big  size  strong-scattering  irregularities 
directly.  However  It  is  sufficient  to  use  the  small-angle 
approximation  according  to  the  big  size  of  irregularities  or  their 
detailes.  In  what  follows  we  consider  the  methods  of  ISP  solutions 
only,  based  on  the  asymptotic  approximation  for  the  forward 
small-angle  scattering.  It  should  be  noted  it  is  particularly 
interesting  for  tomographic  methods  to  consider  such  cases  of 
small-angle  scattering,  when  a  wavelength  of  sounding  radiation  is 
much  less  compared  with  the  scales  of  object  specific  detailes. 
This  is  necessary  condition  to  reconstruct  the  complicated 
internal  structure  of  the  object. 

In  the  beginning  it  is  useful  to  consider  the  ISP  in  the 
high-frequency  limit  i.e.  the  source  field  freguency  is 
essentially  higher  than  the  critical  frequency  of  the  regular 
ionosphere.  In  this  case  the  refraction  index  is  close  to  unity, 
and  the  influence  of  the  regular  ionosphere  may  be  neglected. 
Then,  instead  of  (1),  the  following  equation  should  be  used: 

AE  +  k2E  -  q(r,w)E  =  6(r -rQ)  (4) 

This  differential  equation  is  equivalent  to  the  Lippmann-Schwinger 
integral  equation 

E(r)  =  E0(r)  +  J  G(r-r')  E(r')  q(r',U))  d3r' 


(5) 


where  G(r)=-(4icr)  1exp(ikr)  is  the  Green's  function  for  a  vacuum 
and  E0(r)=G(r-rQ)  is  the  source  field. 

It  is  necessary  to  bear  in  mind  that  the  potential  q  depends 
on  N  and  v  as  follows  (the  usual  condition  for  ionospheric 
radioprobing  v  «  w  is  taken  into  account)  [13 

q  =  4tct  N(r)/(1  +  ir(r)/u)  *  4%r _N(r)(1  -  iv(r)/(i))  (6) 

T.  0  td 


For  the  case  of  small-angle  scattering  the  transmitter,  a 
scattering  object  and  the  receiver  are  situated  along  a  straight 
line  approximately.  Assuming  that  the  direction  of  this  straight 
line  is  closed  to  the  direction  of  z-axis  one  can  obtain  an 
approximate  equation,  describing  the  forward  small-angle 
scattering,  instead  of  (5)  El ,73 

U(r)  =  1  +  |  F(r,r* )  U(r')  q(r',o))  d3r'  =  1  +  FUq  (7) 

Here  U=E/E0.is  the  normalized  field.  Under  the  derivation  (7)  the 
Fresnel's  approximation  was  used  for  the  Green  function  G  either 
as  for  the. sounding  wave  field.  After  some  transformations  we  come 
from  (5)  to  the  formula  (7),  where  the  kernel  F  of  the  equation 
(7)  contains  coordinates  of  the  transmitter  zQ,  the  finite 
scatterer  zs  and  the  receiver  zR 


1  k 

F(r,r')  - - exp  (i~of 

4tcC  ^ 


(p'-  s)2). 


(8) 


®  =  Pr<Z0  _zs^zO  “  ZR^  +  Po^Zs  ZR^Z0  SR^  * 


C  =  (zq  -  2s)(zs  -  zR)/(z0  -  Zr)- 

In  general  case  the  formulae  for  s  and  £  3iave  to  include  the 
integration  variable  z'  instead  of  zs.  But  if  the  object  scale  is 
less  than  longitudinal  Fresnel  resolution  the  varriable  z'  can  be 


changed  on  the  constant  zQ. 

There  is  another  form  of  the  equation  (7),  using  the  complex 
phase  ©  =  In  E,  $  -  ©Q  =  In  U 

©  -  ©0  =  In  (1  +  Fq  exp (©  -  ©Q))  (9) 

which  includes  the  RA  (the  first  iteration  of  the  non-linear 
equation  (9)) 

©  -  ©0  =  In  (1  +  Pq)  «  Pq  or  U  =  exp(Pq)  (10) 

Taking  into  account  (7,10)  and  making  use  of  the  Fresnel's 
transformation  we  have  the  algorithm  for  reconstruction  of  2D 
cross-section  of  potential  q  [1,23 

qz(p,  w)  =  Jq(r,03)  dz  =  (k/2TCC)2  J  V(s)  exp(-ik(s-p)2/2C )d2s. 

(11) 

The  integral  projection  qs(p)  can  be  obtained  after  transformation 
of  field  data.  Por  case  of  the  BA  we  have  V=-4tuC(E-E0)/E0  and  for 
the  case  of  the  RA  V  =  -4tc(;(©-©0)  [1,23.  Clianging  the  direction 
of  the  sounding  wave  incident  ("turning"  the  z-axis)  one  can 
obtain  a  set  of  2D-projections,  which  makes  it  possible  the 
tomography  reconstruction  of  the  3D  structure. 

Determination  of  the  linear  potentials  q(p,w)  (11)  from 
complex  potential  q(r,io)  in  the  case  of  forward  scattering  reduces 
the  ISP  to  a  problem  of  tomographic  reconstruction  -  the  problem 
of  reconstruction  an  object  from  projections.  The  intensive 
development  in  the  past  decade  of  tomographic  methods  for  studying 
the  structure  of  objects  is  to  a  significant  extent  the  result  of 
advances  in  x-ray  tomograpiiy.  The  reconstruction  algorithms  used 
in  practical  x-ray  tomography  are  based  on  a  linear  approximation 
of  the  ray  trajectories.  In  mathematical  respects,  such  problems 
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reduce  to  reconstruction  of  the  attenuation  function  or  refractive 
index  from  group  of  linear  integrals  and  to  reconstruction  of  an 
object  from  its  small-dimensional  projections.  In  addition  to 
x-rays,  practically  all  known  types  of  emissions  and  waves  are  now 
used  for  tomographic  purposes.  The  linear  approximation  frequently 
does  not  provide  good  results  in  tomographic  investigations  using 
optical  and  ultrasonic  waves,  microwaves  and  other  waves. 
Therefore,  methods  for  reconstruction  with  consideration  of 
reflection  and  diffraction  effects  have  been  developed  intensively 
in  recent  years,  and  a  special  term  -  diffraction  tomography  -  has 
appeared  C  8-1 2 , 4 3 . 

The  solutions  of  the  IP  of  reconstruction  of  the  structure  of 
a  scatterer  with  diffraction  considered  in  the  case  of  "forward” 
scattering  obtained  above  can  also  be  applied  to  the  field  of 
diffraction  tomography.  ISP  reduces  to  the  tomographic  problem  of 
reconstructing  a  3D  object  from  2D  projections  qz(p,w)  (11).  By 
rotating  the  object  or  rotating  the  transmitting  and  receiving 
system  in  relation  to  the  object  (axis  s).  It  is  possible  to 
obtain  a  series  of  linear  integrals  of  qz(p,w),  which  can  be  used 
for  tomographic  reconstruction  of  the  refractive  index.  However, 
the  problem  of  reconstructing  the  2D  structure  from  a  set  of 
ID  "cross  sections",  i.e.,  functions  measured  by  one  receiver  as 
the  transmitter  on  a  satellite  moves,  is  also  a  tomographic 
problem.  Furthermore,  the  physical  meaning  of  the  quantities 
measured  by  one  receiver  is  of  integrals  over  the  cross  section  of 
the  spectrum  of  the  irregularity.  A  set  of  such  integrals  makes  it 
possible  to  reconstruct  the  2D  cross  section  of  spectrum,  which  is 
equivalent  to  reconstruction  of  2D  integrated  projection  of 
qz(p,cj).  The  set  of  2D  projections  makes  it  possible  to  perform 
tomographic  reconstruction  of  the  3D  structure. 


The  influence  of  diffraction  reduces  to  the  fact  that  the 
linear  integral  (11)  is  dependent  not  only  on  the  field  "on  the 
ray",  bur  also  on  the  field  in  vicinity  of  the  intersection  point 
of  the  ray  and  recording  plane  (z  =  sR) .  If  the  emission  frequency 
is  high,  then  at  the  limit  w  ->  a>  the  integral  (11 )  should  only  be 
dependent  on  the  field  on  the  ray  (the  eikonal  approximation  with 
a  stright  line  trajectory). 

Let  a  plane  sounding  wave  (pQ=0,  zQ->  oo)  strike  the  object. 
Then,  by  calculating  the  integral  (11)  at  the  of  high  k  (kr^/2C»1 ) 
by  the  saddlepoint  method  (assuming  that  V(s,u))  =  -4^^  in 
accordance  with  smooth  perturbation  method  approximation),  we 
obtain 

qz(p,k)  s  Jq(r,k)dz  -  (k/2-juC)2  V(pR) (2C/k) (-iiu)=  2ik®1(pR). 

(12) 

To  demonstrate,  the  advance  of  the  complex  phase  ©  at  the  object 
in  the  geometrical  optics  approximation  with  a  linear  ray  is  equal 
to  the  integral  of  complex  refractive  index  n(r,w)  =  v  1-q/k2 

©  =  ®0  +  ©  +  ...  =  ik  J  n(r,o))  dz  = 

=  ik  J  v  1-q/k2  dz  -  ik  J  (l-q^k2)  dz 

Prom  this,  for  the  first  approximation  of  the  complex  phase  of  the 
field  we  have  ®1  *  -iqz/2k,  which  agrees  with  (12),  i.e.,  at  the 
high-frequency  limit  linear  integrals  of  type  (11)  are  only 
dependent  on  the  field  on  the  ray.  This  means  that  the 
approximation  used  here  also  contains  the  geometrical  optics 
approximation  with  linear  trajections  at  the  high-frequency  limit. 
The  result  of  the  limitation  on  s  of  the  region  where  the  field  is 
recorded  is  that  the  function  q  (p,w)  reconstructed  from  limited 


data  on  ®1 (s)  will  be  smoothed.  The  null  space  of  the  equation  is 
determined  by  the  parameters  of  the  recording  system  and  defines 
the  Rayleigh  resolution  limit  [13;  then  qz(p,w)  is  an  integral 
over  a  ray  of  finite  "thiclmess". 

To  reconstruct  the  2D  structure  of  ionospheric  irregularities 
according  to  the  transformation  (11 )  we  need  the  information  about 
the  approximate  coordinate  of  the  scatterer  zg.  Let  us  discuss 
briefly  a  method  for  determination  of  a  distance  to  the  scatterer, 
which  is  essential  for  irregularities  reconstruction  by  integral 
transformations.  If  only  an  approximate  value  of  the  distance  is 
known  then  there  appear  errors  in  s  and  C,  and.  hence,  the 
reconstructed  function  is  destorted.  Let  us  consider  an  example 
where  a  source  is  placed  on  board  of  a  satellite  and  the 
corresponding  distortions  arising  in  this  case.  The  distance 
between  the  satellite  and  the  receivers  (zQ  -  zR)  is  known  quite 
precisely,  but  in  determination  the  distances  "source  -  scatterer" 
(Zq  -  z  )  and  "scatterer  -  receivers"  (zs  -  zR)  a  systematic  error 
A  is  present:  zg  =  zg  -  A.  It  may  be  shown  that  reconstruction 
errors  are  small  for  such  values  of  A,  which  correspond  to  the 
longitudinal  resolution  of  the  measuring  system.  This  conclusion 
is  also  confirmed  by  numerical  simulations  [13. 

It  easy  to  show  that  when  the  error  A  is  less  than  the  longi¬ 
tudinal  resolution,  A  «X  2C2/(max  s)2,  the  reconstructed 
function  is 

qz(p)  “  qz(2/V  y/v  exp(12Ciy^  (y2_  x2))’  (13) 

where  ux  =  1  -  A/(zg-zR),  =  1  +  A/(zQ-zs).  This  equation  makes 
it  possible  to  determine  the  distance  to  the  scatterer  by 
performing  reconstruction  with  different  values  of  A.  At  the  true 
location  of  the  scatterer  (A=0)  the  "phase  front"  curvature  of  the 


function  under  reconstruction  changes  to  the  opposite  with 
respect  to  each  coordinate.  These  problems  are  considered  in 
detail  in  [13  and  partly  in.  [63-  There  one  can  find  the 
generalisation  of  the  method  for  the  ionosphere  with  the 
stratified  background;  the  method  for  determination  of 
inhomogeneities  parameters  ("mass",  space  coordinates,  size  and 
other  moments  of  potential  q)  by  means  of  a  small  number  of 
receivers  and  solutions  of  IP  due  to  data  of  back-  scattering. 

It  is  possible  to  consider  the  more  general  case,  namely,  the 
strong  scattering  case.  For  this  case  in  accordance  with  (7,10)  we 
reconstruct  the  projection  of  the  product  of  the  potential  and 
the  normalized  field  Q=qU,  namely, 

J  qU  dz  =  Qz(p) 

It  should  be  stressed  that  in  general  the  normalized  field  UCq3 
depends  on  the  potential  and  on  the  direction  of  the  sounding  wave 
incident.  So  it  is  impossible  to  solve  the  problem  of  the  product 
qU  reconstruction  by  means  of  the  direct  transformation  of  a  set 
of  projections  (14).  It  is  possible  to  suggest  Iteration 
procedures  of  the  ISP  solution  to  the  equation  (14)  directly.  But, 
if  the  scale  of  ionospheric  inhomogeneities  more  than  sounding 
wavelength  in  103-  104  times,  it  is  even  complicated  to  solve  the 
direct  scattering  problem.  So  for  such  cases  It  is  necessary  to 
use  asymptotic  methods.  The  special  asymptotic  methods  of  the 
solution  of  the  ISP  for  the  case  of  strong  and  large  scatterers 
were  developed  in  [1,183.  For  example  the  following  asymptotic 
representation  of  normalized  field  was  obtained 

i  3  if  d2q(p+z)  _o 

U(p,z)  =  exp[-  -  fq(p+z')dz  +  — p  (z-z') — — y-  dz'  +  0(k  -)3 

2k  4k2  dp- 


(15) 


This  formula  is  correct  for  infinitly  smooth  potential.  The 
corresponding  additional  terms  in  the  sum  (15)  will  appear  in  the 
presence  of  derivatives  discontinuities.  The  formula  (15)  permits 
us  to  construct  the  simple  iteration  procedure  of  the  ISP 
solution,  when  it  is  not  necessary  to  solve  numerically  neither 
the  direct  scattering  problem,  or  the  ISP.  The  given  first  guess 
of  the  potential  q1  (a  priori  information  or  transformation  under 
the  condition  of  weak-scattering)  defines  the  approximation  of  the 
normalised  field  functional  Utq.,]  according  to  (15).  Using  this 
approach  it  is  possible  to  solve  the  usual  tomographic  problem  of 
the  reconstruction  q(r)  (14)  with  known  "weight"  U  due  to  Qz(p) 
under  the  different  directions  of  the  sounding  wave  incident, 
later  on  the  iteration  procedure  with  the  potential  approximations 
obtained  is  repeated. 

The  formula  representation  (15)  permits  us  to  formulate  the 
conditions  of  applicatability  of  the  weak-scattering 
approximations.  According  to  this  condition  the  first  term  ("k-1 ) 
and  all  successive  terms  of  the  expansion  (15)  have  to  be  small. 
The  restriction  on  the  first  term  gives  us  the  well-known  • 
condition  of  the  BA  applicatability. 

<Fm/2k  «  1  (16) 

Here  q  is  the  typical  (mean)  value  of  the  potential,  rm  Is  the 
maximum  size  of  the  inhomogeneity.  In  the  case  when  the  wavelength 
of  the  sounding  radiation  X=2m  and  the  value  of  disturbances  of 
the  Ionosphere  electron  density  N  ~  1 Q1 1  el/m3  the  condition  (16) 
leads  to  the  inequality  rm«  2  km.  The  disturbances  of  the 
electron  density,  scattering  radiowaves,  have  to  be  measured 
relative  to  the  background  of  the  regular  stratified  ionosphere. 
So,  if  the  background  value  N  1 01 1  el/m3  and  the  disturbances  are 


sufficiently  strong  (10%),  the  BA  is  valid,  up  to  the  sizes  of 
dozens  km.  But  for  the  case  of  the  main  ionosphere  maximum  (N  ™ 
lO1^  nf3)  the  scale  restriction  is  significant  and  it  is  necessary 
to  take  into  consideration  the  first  term  of  (15).  The  second  term 
of  (15)  for  simple  scatterers  without  internal  structure  is  less 
^  ^  than  the  first  one.  In  the  case  when  there  are  some 
internal  details  inside  the  scatterer  with  the  scale  a»A  and  the 
typical  value  of  the  potential  q',  the  second  term  can  be  more 
than  the  first.  Hence  the  other  necessary  condition  of  the  BA  is 

the  condition  that  the  second  term  with  the  transverce  derivative 
is  small 

q'rm/(4k2a)  «  i.  (17) 

In  many  practical  cases  It  is  sufficient  to  take  into  account 
the  first  term  of  (15)  only,  but  the  scatterer  is  not  the  Bom's 
scatterer  according  to  (16)  and  the  inequality  (17)  is  valid.  Such 
conditions  allows  us  to  get  the  simple  analytic  formula  El  ,143, 
which  connects  the  function  Qs  reconstructed  due  to  the 

experimental  data  and  the  projection  of  the  scattering  potential 
q3  (15). 

Qs(p)  =  2ik  [exp(-Iqs(p)/2k)  -  1]  (18) 

This  formula  is  the  basis  of  the  tomography  reconstruction  for  the 
most  of  ionospheric  applications.  Really,  the  inhomogeneities  with 
the  scale  about  ten  or  more  kilometers,  the  deviations  of  which 
are  about  10%  with  respect  to  the  background,  are  the  strong 
inhomogeneities  (the  left  hand  side  of  (16)  Is  mors  or 
approximately  equal  to  6,  when  k=  2  m),  but  the  inequality  (17)  is 
valid  for  the  internal  details  with  the  scale  more  than  hundreds 
of  meters  and  with  the  electron  density  disturbances  about  100% 
compared  with  the  background. 


The  uniqueness  of  solutions  of  the  similar  ISP  was  considered 
many  times  earlier  [1-31.  Some  exact  results  are  known.  For  the 
problems  in  question  the  uniqueness  of  the  transformation  qz  (11) 
and  the  following  reconstruction  (14,18)  of  projections  can  be 
easily  proved  for  the  finite  volume. 

The  results  described  above  allow  us  to  get  2D  projections  of 
localised  ionospheric  inhomogeneities  with  the  help  of  the  moving 
satellite  transmitter  and  the  transverse  arrow  of  Earth  receivers. 
If  there  are  some  arrows  of  receivers  separated  by  distances  about 
hundreds  kilometers  along  the  line  of  the  satellite  flight,  it  is 
possible  to  obtain  some  2D  projections.  Having  a  set  of  2D 
projections  one  can  reconstruct  the  3D  structure  also.  Such 
experiment  (with  some  arrows  of  receivers)  is  rather  complicated 
and  it  was  not  realized  yet.  Nowadays  experiments  connected  with 
the  reconstruction  of  the  2D  localized  inhomogeneities  were 
carried  out  [1,53.  This  is  the  tomography  problem  also  because  it 
is  possible  to  reconstruct  the  2D  inhomogeneity  projection  or  the 
2D  cross-section  of  its  spectrum  after  the  corresponding 
transformations  by  means  of  a  set  of  ID  field  records  in  each 
receiver. 

We  will  consider  the  examples  of  numerical  calculations  by 
means  of  eleborated  programs. 

Using  the  first  numerical  example,  we  will  illustrate  the 
necessity  for  considering  diffraction  effects  in  RT  reconstruction 
of  objects  with  dimensions  rm  comparable  to  the  dimensions  of  the 
Fresnel  zone  vT£.  As  was  demonstrated  above  (12),  diffraction 
effects  are  insignificant  where  kr^/2C  »  1  and  can  significantly 
distort  the  results  where  krj^/2C  <  1.  Fig.  2a  shows  qz(p)» 
consisting  of  two  irregularities  -  two  Gaussians  (~  exp(-p2/r|) 
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with,  parameters  (r^/v  A,C  =  5,  r2/v  XC,  -  1/3,  for  the  A-2m, 
C=200km,  2r1  =  6.3km,  2r2  =  420m).  If  we  reconstruct  the  complex 
phase  after  the  emission  passes  through  one  large  irregularity, 
then  according  to  (12)  Im  fc,  is  proportional  with  a  high  accuracy 
to  the  function  qz(p)  from  a  large  irregularity.  However,  the 
imaginary  part  of  the  complex  phase  ®1  (Fig.  2b)  from  a  pair  of 
irregularities  bears  little  resemblance  to  Fig.  2a  due  to 
diffraction  effects  in  small  irregularity.  Fig.  3  shows  the 
similar  picture  for  the  two  irregularities  with  parameters 
(r/vTC  =  5,  r2/v~XC  =  1/8).  The  pair  of  irregularities  is 
reconstructed  very  accurately  in  the  case  where  the  Fresnel 
inversion  (11)  is  used. 

For  reconstruction  of  the  structure  of  ionospheric 
irregularities  the  questions  about  sampling  data  and 
reconstruction  procedure  are  significant.  The  sounding  radio 
signals  are  always  recorded  by  point  receivers.  The  signal  from 
each  receiver  recording  the  emission  of  the  moving  transmitter  on 
the  satellite  is  also  sampled  by  an  electronic  circuit  with  a 
fixed  interrogation  frequency.  Moreover,  the  integral 
transformations  as  a  result  of  which  the  structure  of  the 
scattering  irregularities  is  reconstructed  can  be  realized 
numerically  on  a  computer  only  after  the  appropriate  digitization. 
Therefore,  both  the  form  of  the  data  recorded  and  the  nature  of 
the  numerical  processing  and  the  reconstruction  of  the  objects 
lead  to  the  necessity  for  performing  purely  discrete  operations. 
However,  in  theoretical  examination  of  the  solution  of  the  inverse 
scattering  problem  the  immediate  switch  to  discrete  formulas  is 
hardly  advisable  in  the  majority  of  cases.  Usually,  the  discrete 
analogues  to  the  reconstruction  formulas  are  significantly  more 


cumbersome,  and  besides  the  discrete  form  hinders  analysis  of  the 
results  obtained.  On  the  other  hand,  the  switch  from  continual 
relations  to  their  discrete  analogues  does  not  cause  difficultes 
in  the  majority  of  cases  and  is  performed  in  a  unique  manner.  So, 
continual  formulas  can  be  considered  a  compact  notation  of  their 
discrete  analogues  for  further  practical  reconstruction. 

Let  us  consider  the  switch  to  discrete  analysis  of  the 
reconstruction  procedure  and  the  questions  associated  with  this 
switch  using  the  transformations  (8),  (11).  It  should  be  noted 
that,  the  methods  for  digital  processing  of  signals  and  fields 
have  been  developed  rather  well  and  covered  in  the  literature 
[19,20]. 

First,  we  switch  in  (8,  (11)  to  the  dimensionless  variables  P 
=  p/vTC  ;  S  =  s/vTC  ,  normalized  to  the  radius  of  the  Fresnal 
zone.  Intriducing  the  quantities 
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from  (11)  we  obtain  the  pair  of  integral  Fourier  transforms  for 
them  [1  ] 


F  (P)  =  X  FV(S)  e  2ilcSP  d2S 

T, 

F  (S)  =  I  F  (?)  e+2l7tSP  d2p 

V  SI 


(20) 


By  virtue  of  the  definition  Fq(P)  is  a  finite  function  due  to 
the  finiteness  of  qz(p).  Let  the  dimensions  of  the  object  not 

exceed  the  limits  of  the  segment  [-Pqx,  Pqx)  ^  x  ^-P0y’ 

P  ]  in  axis  y,  i.e.,  the  carrier  of  the  function  F  (P)  is 
Oy 

contained  in  a  fixed  rectangle.  Then,  FV(S)  will  be  an  infinite 


and  analytical  function  with  a  finite  spectrum  F  (P).  According  to 
the  Kotel’nikov  theorem  FV(S)  can  be  represented  in  the  form  of  an 
infinite  series  in  terms  of  sample  functions  with  the  sampling 
Intervals  ASX  =  P^/2,  ASy  -  P0y/2:mm 

VW  =  sine (2®^-®^)  slneCPtfpySy-imiy) 


Substituting  this  representation  in  (20)  and  integrating  with 

respect  to  S  ,  S  ,  we  obtain  the  relation 

y  (21 ) 

“y  -2iic[P  AS jL  +  P^s jL ]  1 

F„(P)  =  F  ( - ,  - )  e  x  x  x  3r  y  y  /?p  ?p  v-1 


Fv( - ,  - )  e 

®x  "V  2P0x  2P0y 


^Ox  2^0y^. 


Using  the  concept  of  signals  which  are  really 

indistinguishable  at  the  level  e  [213.  it  is  advisable  to 

introduce  a  finite  function  which  is  really  indistinguishable  at  a 

-»  w 

fixed  level  from  FV(S).  Let  the  carrier  of  this  function  be 

contained  in  the  rectangle  [-S0x,S0x]  x  [-S0y,S0y].  Then, 
according  to  the  Landau-Pollack  theorem  [21,22]  the  approximate 

dimension  of  the  set  of  all  functions  FV(S)  finite  at  a  fixed 

level  with  a  carrier  in  the  rectangle  [-S0z,S0x]  x  [-SQy,S0y]  and 
with  a  finite  spectrum  supp  Fq(P)  €  C-P0x,  PQx]  *  [-PQy,  PQy]  is 
close  to  =  (2S0x><2P0x)  x  (2S0yx2P0y).  In  other  words,  due 
to  the  finite  accuracy  in  measurement  of  the  signals  it  can 
assumed  that  both  the  signal  and  its  spectrum  are  finite.  The 

limitation  on  the  carrier  of  the  spectrum  leads  to  a  finite 
sampling  frequency  and 

Y?),  "x  =  1/2S0x-  ATy  =  '/2S0y. 
here,  the  number  or  sampling  Intervals  on  the  Cartesian  axes  Is 
equal  to  Nx  =  =  =  2P0x/APz.  Ny  =  2Poy/APy. 

From  this,  assuming  Px  =  nxAPz,  py  =  n  AP  ,  It  Is  simple  to  switch 


from  (19),  (20)  to  a  pair  of  discrete 


F/CAP^APy)  H  Fq 


=  F 


q 


transformations  for  Fy  = 


1  Nx/2-1  Ny/2-1 

W  mx=-Nx/2  my__Ny/2 


V^iy 


-2iiu( 


Vi  +  Vy 


1L 


N 


(22) 


VW 


Nx/2-1  Ny/2-1 

Vnx/2  V"Ny/2 


^q^x*1^ 


mxnx+W 


Later  on,  the  fast  Fourier  transformation  algorithms  are  used  in 
the  numerical  modeling  and  the  numbering  of  the  sums  over  n^,  n^ 
is  shifted  from  0  to  Nx  -  1 ,  and  likewise  for  ny 

The  dimensions  of  the  reception  region  SQx,  SQy  and, 
consequently,  the  sampling  intervals  APX,  APy  also  determine  the 
resolutions.  To  demonstrate,  the  minimum  resolvable  interval  Sx  in 
x  is  equal  to  5X=  7~X C,  APX=  YTC/2S0x  =  =  A,C/2Sm.  A  limitation  of 
the  reception  region  always  limits  the  resolution  or  the  sampling 
interval  and  the  function  with  a  finite  spectrum.  Therefore, 
generally  speaking,  with  a  limited  reception  region  it  is  possible 
to  switch  immediately  from  (21)  to  the  finite,  discrete  Fourier 
transform  (22). 

The  switch  from  other  continual  reconstruction  formulas 
containing  the  Fourier  and  Fresnel  transforms  to  their  discrete 
forms  is  also  made  according  to  this  scheme. 

The  practical  realization  of  reconstruction  of  the  structure 
is  governed  to  a  great  extent  by  the  resistance  of  the 
reconstruction  procedure  to  the  various  types  of  noises  and 
distortions  which  unavoidably  appear  in  measurements  and 
processing  of  experimental  data.  Let  us  illustrate  the  influence 


of  noises  and  distortions  on  reconstruction  of  the  structure  on 
the  examples  of  numerical  simulation.  Briefly  discribe  a  modeling 
scheme  using  the  example  of  two-dimensional  Fresnel 
reconstruction. 

The  modeling  scheme 

qz(p)  =>  Fq(P)  ->  PV(S)  =>  V(S)  =>  .....  =>  (^(p) 

! 

j 

noises,  distortions 

The  function  qz(p)  characterizing  the  two-dimensional  structure  of 
the  electron  concentration  irregularities  and  the  effective 
collision  frequency  is  defined.  Then,  the  quantity  Fq(P)  is 
calculated  using  formula  (19).  After  performing  the  discrete 
Fourier  transform  (22)  corresponding  to  (20),  Fy(S)  is  obtained 
and,  thus,  the  field  V(S)  in  a  discrete  grid.  This  stage  of 
modeling  is  related  to  the  direct  scattering  problem,  where  it  is 
assumed  that  the  scattering  is  calculated  well  within  the 
framework  of  the  corresponding  approximations.  Further,  the  noise 
can  be  added  to  the  data  on  the  field  obtained,  after  which  the 
inverse  discrete  Fourier  transform  (22)  corresponding  (20)  is 
performed,  and  the  result  of  reconstruction  qz(p)  is  obtained  as  a 
result.  The  influence  of  distortions  is  modeled  in  asimilar 
manner:  in  the  inverse  transformation  stage,  the  parameters  of  the 
transformation  are  changed  or  ’’distorted”. 

The  results  of  numerical  modeling  of  the  influence  of  noises 
have  demonstrated  the  stability  of  the  reconstruction  procedure. 
Fig. 4  shows  two  actual  irregularities  of  qz(p)  -  two  Gausslans 
with  various  dimensions  on  a  64x64  grid.  The  total  size  of  the 
image  frame  in  units  of  the  scale  of  the  Fresnel  zone  VTC  is 
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equal  5x5.  After  calculating  the  complex  phase  of  the  field 
scattered  by  these  irregularities,  the  data  on  V  were  disturbed  by 
noise.  Then,  reconstruction  was  performed;  the  quality  of  it  was 
sufficient  even  with  comparatively  large  errors  in  the  measured 
field  [23].  For  example.  Fig. 5  shows  the  result  of  reconstruction 
of  these  irregularities  according  to  data  on  with  additive 
complex  Gaussian  noise  having  a  variance  of  0.05  of  maximum 
amplitude  of  the  change  in  the  real  and  imaginary  parts  of  complex 
phase  of  the  field.  Fig.6.  illustrates  reconstruction  with  double 
the  noise  (a  variance  level  of  0.1).  Another  numerical  example: 
actual  irregularities  -  "cylinder”  and  "parabola"  are  represented 
in  Fig.T.  The  result  of  reconstruction  of  this  model  structure 
with  additive  noise  (a  variance  level  of  0.05)  is  shown  in  Fig.8. 
Fig.9  illustrates  the  reconstruction  with  double  noise  (0.1).  If 
we  perform  preliminary  processing  of  the  data,  then  even  higher 
noise  levels  will  not  noticeably  influence  the  results  of 
reconstructions. 

Numerical  estimation  of  the  level  of  influence  of  the  noises 
is  of  interest  for  practical  applications.  It  is  advisable  to 
characterize  the  level  of  influence  of  the  noises  in  the  metric  G 
and  L2.  In  a  discrete  reconstruction  procedure,  it  is  better  to 
use  their  discrete,  normalized  analogues,  i.e.  we  will  estimate 
the  difference  between  the  function  f  reconstructed  from  noisy 
data  and  real  function  f  by  the  numbers 


P0(f.?> 


rv 


max|f.-f.| 

i 

max  |  ^  | 


P  2(f,f) 

if- 


i  1  1  1/2 

-  ] 

2  f? 


(23) 


where  summation  and  selection  of  maximum  are  performed  for  all 
increments  of  recorded  function. 


The  dependence  of  the  normalized  disturbance  of  the 
reconstructed  two-dimensional  structures  (i.e.,  the  difference 

rv 

between  the  true  real  qz(p)  and  reconstructed  qz  in  the  metrics 
p  and  p  0)  on  the  amplitude  A  of  normalized  (to  the  maximum 

C  jd 

amplitude  of  the  change  in  the  field)  complex  Gaussian  noise  is 
shown  in  Tabl.1  for  the  irregularities  in  Fig.?.  It  is  clear  that 
the  normalized  deviations  of  the  reconstructed  two-dimensional 
structures  are  comparable  to  the  noise  level.  The  absolute 
disturbances  of  the  real  part,  although  the  true  function  is  real. 

Distortions  of  irregularities  reconstructed  by  radio 
tomographic  methods  appear  due  to  imprecise  determination  of  the 
distance  to  the  irregularities  and,  as  a  result,  incorrect 
knowledge  of  the  given  parameters  of  the  integral  transformations. 
The  influence  of  distortions  was  modeled  by  the  scheme  described 
above.  The  model  in  Fig.  10  (two  irregularities:  "cos”  and  "cosA2" 
ellipse)  is  the  object  of  reconstruction.  The  imprecisely 

determined  coordinate  of  the  scatterer  zg  =  zg  -  A  with  an  error  A 
was  substituted  in  the  parameters  of  the  Fresnel  transformation. 
Numerical  values  zQ  -  zg  =  700  km,  zg  -  zR  =  300  km.  A,  =  2  m 
typical  for  satellite  sounding  were  used  in  the  calculations.  The 
dimension  of  the  Fresnel  zone  here  V  AC  *  0,65  km,  the  frame 

dimension  was  6, 4x6, 4  km  and  the  transverse  resolution  6  was  of 
the  order  of  a  hundred  meters.  Fig.  11 -15  provide  reconstructed 
two-dimensional  structures  with  various  values  of  the  error  in 
determination  of  the  coordinates  A  measured  in  kilometers;  A  is 
successively  equal  to  -10,  +20,  -20,  +30,  +40  km.  It  is  clear  that 
an  error  of  10  km  has  little  influence  on  the  reconstruction. 

Errors  A  of  the  order  of  20  km  distort  the  object,  but  the 

characteristic  details  of  the  object  are  reconstructed  fairly. 
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Only  errors  exceeding  30-40  km  noticeably  distort  the 

reconstructed  irregularity  and  may  lead  to  Incorrect  conclusions  on 
the  structure  of  the  object.  The  values  of  the  errors  In 
reconstruction  in  normalized  metrics  are  shown  In  Tabl.2.  The 
modeling  was  performed  for  two  Gaussians.  The  analogy 

dependence  of  error  reconstruction  on  the  A  is  shown  in  the  graphs 
in  [1],  The  reconstruction  errors  have  a  clearly  pronounced 
minimum  and  are  minor  (<  0,1  +  0,2)  within  ±  (5  -  10)km.  We  will 
note  that  dimension  of  this  minimum  agrees  in  order  of  magnitude 
with  the  longitudinal  resolution  5Z  of  the  system;  here,  the 
transverse  resolution  6  ~  X/Q  ~  100  m,  the  aperture  angle  6  ~0,02 
h  6  ~  A/62  ~  5  km.  Therefore,  determination  of  the  distance  to 

scatterer  with  accuracy  of  the  order  of  the  longitudinal 
resolution  of  the  recording  system  is  sufficient  for  high-quality 
reconstruction.  The  method  for  determining  the  distance  according 
to  the  change  in  the  curvature  of  the  "phase"  of  the  reconstructed 
function  makes  it  possible  to  achieve  this  accuracy  in 
insignificant  noises  [1]. 

On  reconstruction  of  qz  using  the  formula  (18)  the  main 
obstacle  is  that  Qqz  is  a  periodical  function  of  qz,  therefore,  a 
real  part  of  qz  is  not  unambigiously  defined  by  Qz,  which  we  get 
from  experimental  data.  Since  that,  making  use  the  formula  (18) 
one  cannot  immediately  reconstruct  scatterers  that  exceed  47uk.  It 
should  be  noted  here  that  it’s  much  more  than  BA-based 
reconstruction  procedures  are  applicable  to.  Namely,  we  have  found 
using  numerical  simulations  that  BA  holds  for  scatterers 
satisfying  the  inequality 


qz/2k  <0.1. 


In  principle,  one  can  reconstruct  more  intensive  scatterers, 
assuming  that  Qz  is  zero  if  qz  crosses  ick  level,  or  if  any 
additional  information  about  the  scaterer  structure  is  available. 

We  have  derived  an  algorithmic  procedure  for  finite 
scatterers,  which  allows  one  to  reconstruct  scatterers  exceeding 
4iuk.  We  have  seen  that  the  main  obstacle  which  restricts  this 
procedure  capabilities  is  the  finite  number  of  receivers  which 
register  the  scattered  field.  From  the  formal  mathematical  point 
of  view,  that  simply  means  we  get  only  finite  number  of  Fourier 
expansion  terms  ,  i.e.  we  deal  with  the  truncated  Fourier  series 
of  Qz. 

Thus  to  obtain  q„  we  have  to  derive  it  from  a  truncated 
Fourier  expansion  of  Qz,  which  in  turn  is  a  nonlinear  transcendent 
function  of  qz.  This  makes  us  use  the  numerical  simulations  to 
investigate  the  domain  of  applicability  of  this  procedure,  because 
immediate  analysis  of  that  algorithm  is  extremely  difficult. 

We  have  performed  such  an  investigation  of  this  procedure 
applying  it  to  various  scatterers.  The  examples  of  such  numerical 
simulations  of  reconstruction  of  ID  projections  of  2D  scatterers 
are  presented  on  Fig. 16-1 8.  The  scatterer  -  receiver  distance  was 
assumed  to  be  300  km,  the  sounding  wavelength  was  2m.  Thick  curves 
correspond  to  initial  projections  qz  which  had  to  be 
reconstructed,  thin  curves  are  the  reconstructed  images  of  qz 
using  the  procedure  described  above.  Just  for  checking,  the  Qz 
patterns  are  also  shown  (small  amplitude  fast  oscillating  curves 
of  normal  thickness).  On  the  Figure  16,  the  scatterer  is  not  too 
strong,  so  two  curves  simply  coincide  (  the  projection  is  exactly 
reconstructed  ).  On  the  Fig. IT-1 8,  where  stronger  patterns  are 
given,  one  can  see  the  discrepancy  between  exact  image  and 


reconstruction.  One  can  see  very  easy  that  mistakes  in  the 
procedure  work  arise  mostly  on  those  fragments  of  scatterers  where 
the  transverse  derivative  !dqz/dxi  =  !  dQz/dx !  is  big  enough,  so 
the  oscillations  of  Q„  become  so  fast  that  CL  is  no  longer 
approximated  by  the  truncated  Fourier  expansion  well  enough.  They 
may  cause  the  loss  of  zeros  of  Qz,  which  correspond  to  increments 
or  decrements  of  4ick.  The  numerical  simulations  using  wide  variety 
of  scatterers  have  confirmed  this  conclusion.  Moreover,  they  have 
shown  that  only  the  transverse  derivative  absolute  value  idqz/dxi 
is  critical  for  the  correct  reconstruction.  Of  course,  it  causes 
the  limitation  on  the  maximal  intensity  of  the  scatterer,  since 
for  continuous  function  q  (x)  ! qz i  <  a  idqz/dxi.  For  this  reason, 
the  problem  to  calculate  this  boundary  value  for  idqz/dxi  had  also 
been  posed.  The  value  which  we  get  estimating  the  supremum  of  the 
truncated  Fourier  expansion  derivative  (it  is  bounded  as  a 
derivative  of  polynom)  is  much  more  than  actual  value.  Fig.  19 
shows  the  normalized  derivative  value  times  plotted  versus  N, 
where  N  is  the  number  of  terms  retained  in  the  Fourier  series.  We 
assumed  that  well-approximating  reconstruction  it  that  one  which 
has  the  maximal  normalized  discrepancy  between  exact  and 
reconstructed  images  less  than  0.1. 


The  description  of  Fortran  Programs 


System  Requirements 


-  Computer: 

-  Operating  System: 

-  Memory: 


-  Hard  Disk  Space: 

-  Software: 


IBM  AT  or  compatible  (with  coprocessor) 
MS-DOS  or  PC-DOS  version  3.0  and  later 
at  least  Extended  memory  4  Mbytes 
(depends  on  geometry  and  type  of 
approximation  reconstructed  function) 

4  Mbytes 

compiler  1.4e  and  linker  2. 2d 
NDP -FORTRAN-386 (c )  MicroWay  or  later 


1.  Program  <DT_DIR.for> 

This  program  solves  direct  problem,  namely,  determines  complex 
field  V(3)  by  means  of  Fast  Fourier  Transformation  (FFT)  based  on 
the  model  structure  of  irregularities. 

Input  parameters  and  files: 

(the  determination  of  the  irregularities  parameters) 
file  <mod_par . tsk>  contains  the  following  parameters: 

«  line  1  in  file 


’wave  length  (km):’ 

0.002 

(DLV  in  program) 

line  2 

’distance  to  satellite  (km):’ 
line  3 

1000. 

(ZZO  in  program) 

'distance  to  irregularities  (km):* 

300. » 

(Z  in  program) 

These  parameters  you  can  change  in  the  file. 

NY  -  number  of  receivers 

MNY  -  the  degree  of  2  for  NY  (namely,  2AMNY=NY) 

NX  -  number  of  discrete  intervals  on  the  axis  X 
MNX  -  the  degree  of  2  for  NX  (namely,  2AMNX=NX) 

RMN  -  the  radius  of  the  Fresnel  zone. 

PX.PY  -  the  limits  of  segments  [-PX/2,  PX/2]  in  axis  X  and 
[-PY/2,  PY/21  in  axis  Y;  object  (model  structure  not 
exceed  these  limits;  variable  P) 

SX,SY  -  the  limits  of  rectangle  C-SX/2,  SX/2]x[-SY/2,SY/2] 
(variable  3) 
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YR  -  the  size  of  the  receiving  system  (km) 

X02  -  the  size  of  aperture  ’’synthesized"  by  satellte  (km) 

DXX,  DYY  -  the  size  of  frame  (km*km) 

NMOD  -  number  of  model  structure  (1-10): 

1  -  the  gauss Ians 

2  -  ’’homogeneous"  ellipse 

3  -  "parabolic"  ellipse 

4  -  ”  2  +  3  " 

5  -  "  2  "  with  absorption  10% 

6  -  "  3  "  with  absorption  20% 

7  -  "cos"  ellipse 

8  -  "cosA2"  ellipse 

9  -  "  7  +  8  " 

10  -  arbitrary 

Nmod  is  introduced  from  the  screen 
Const5=0.1  (namely, 10%)  -  constant  for  NM0D=5 
Const6=0.2  (20%)  -  constant  for  NM0D=6 

Z1,Z2,Z3,UG  - 

A1,B1,A2,B2,X0,Y0,Y1,A3,FPC,FPS  parameters  for 

AELP.BELP 

APAR.BPAR  models 

AELP1 ,BELP1 ,APAR1 ,BPAR1  - 

RN  -  complex  field 

Subroutine  <SMAXMI>  calculates  min  and  max  of  array 
Subroutine  <FFT>  calculates  1-D  Past  Fourier  Transformation 

Output  files: 

<M0DEL.GRD>  -  contains  the  model  structure  of  irregularities  and 
the  number  of  the  model  (NMOD)  and  constans  VMD  and 
ZRD  (for  determination  of  reconstruction’s  errors) 
<V_real.grd>  -  the  array  of  real  parts  of  the  field  RN 
<V_imag.grd>  -  the  array  of  imaginary  parts  of  the  field  RN 
Execution 
f77  dt_dir.for 
RUN 

dt_dir.exp 


2.  Program  <DT_INV.for> 

This  program  solves  inverse  scattering  problem  for  weak- 
scattering  irregularities  in  the  case  when  noise  is  present  in  the 
measurements  of  the  field  and  calculates  the  errors  of  the 
reconstruction  in  metric  C  and  12  in  dependence  on  level  of  noise. 

Input  parameters  and  files: 

Parameters  <  NX,  MNX,  NY,  MNY,  PX,  PY,  SX,  SY  >are  similar  to  same 
parameters  of  program  <DT_DIR.F0R>. 

PILES  (from  program  <DT_DIR.for>): 

<M0DEL.GRD>  -  input  file  of  model  structure. 

<Y_real.grd>  -  the  array  of  real  parts  of  the  field  RN. 
<V_imag.grd>  -  the  array  of  imaginary  parts  of  the  field  RN. 
Parameter  <Error>  -  the  level  of  noise  (0  <  error  <  1 ) 

Error  is  introduced  from  the  screen. 

ConstS,  Const6  -  are  similar  to  same  parameters  of  program 
<DT_DIR. for> 

Subroutine  <SMAXMI>  calculates  min  and  max  of  array 
Subroutine  <FFT>  calculates  1-D  Fast  Fourier  Transformation 
Function  <RAN>  -  determinates  random  values  in  [0,1]. 

K33  -  the  constant  for  function  <RAN> 

Output  files: 

REC_er.GRD  -  the  reconstruction 

er_rier.dat  -  contains  the  errors  of  reconstruction  in  metric  C 
and  L2 

Same  errors  are  shown  on  the  screen. 

Execution 

f77  dt_inv.for 

RUN 

dt_inv.exp 

3.  Program  <DT_ERR.for> 

This  program  solves  inverse  scattering  problem  for  weak- 
scattering  irregularities  in  the  case  when  distortions  are  present 
in  the  measurements  of  the  field  and  calculates  the  errors  of  the 
reconstruction  in  metric  C  and  L2  in  dependence  on  distortions. 
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Input  parameters  and  files: 

Parameters  <  NX,  MNX,  NY,  MNY,  PX,  PY,  SX,  SY  >  are  similar  to  same 
parameters  of  program  <DT_DIR.F0R>. 

FILES  (from  program  <DT_DIR.for>): 

<M0DEL.GRD>  -  input  file  of  model  structure. 

<V_real.grd>  -  the  array  of  real  parts  of  the  field  RN. 
<V_imag.grd>  -  the  array  of  imaginary  parts  of  the  field  RN. 
<mod_par.  tsk>  is  similar  the  same  file  in  program  <DT_DIR.for>. 
Parameter  <DEL>  -  the  error  of  irregularities  coordinate  in  km 
DEL  is  introduced  from  the  screen. 

Subroutine  <FFT>  calculates  1-D  Fast  Fourier  Transformation 

Output  files: 

REC_cr.GRD  -  the  reconstruction 

er_r!cr.dat  -  contains  the  errors  of  reconstruction  in  metric  C 
and  12 

Same  errors  are  shown  on  the  screen. 

Execution 
177  dt _ err. for 

RUN 

dt_err.exp 


4.  The  solution  of  the  statistical  radiotomography  problem 
(statistical  inverse  scattering  problem) 

Rather  frequently,  the  ionosphere  contains  entire  regions 
filled  with  large  numbers  of  irregularities  of  various  dimensions. 
Such  a  state  is  typical  for  the  equatorial  and  polar  regions, 
especially  at  night.  In  this  case  it  is  advisable  not  to 
reconstruct  individual  realization  of  such  disturbed  ionosphere, 
but  rather  to  pose  the  problem  of  reconstruction  of  the 
statistical  characteristics  of  the  turbulent  random  ionosphere 
such  as  correlation  function  or  spectrum  of  fluctuations.  So,  such 
reconstruction  of  random  ionosphere  statistical  properties  by  the 
measured  field  statistics  it  is  reasonable  to  name  the  statistical 
ISP  or  the  statistical  RT  problem  [24,25]. 

We  consider  the  reconstruction  of  ionosphere  electron  density 
fluctuations  spectrum  by  means  of  satellite  radioprobing  data. 
Integral  equations,  which  relate  the  measured  field  to  the  medium 
structure  are  probably  the  most  convenient  and  adequate 
mathematical  technique  for  the  tomography.  In  this  paper  we  shall 
show  how  to  derive  such  equations.  For  simplicity  it  is  reasonable 
to  use  two  frames  of  reference.  The  first  "global”  frame  of 
reference  r=(x,y,z)=(p,z),  its  origin  is  related  to  the  receiving 
system  as  in  the  previous  section  (Fig.  1 ) .  The  origin  of  the 
second  frame  of  reference  (X,y,Z)= (P,Z)  is  reasonable  to  locate 
at  one  of  the  transmitter  positions  and  the  Z  axis  to  direct  to 
the  center  of  the  receiving  system.  As  the  transmitter  moves  it  is 
convenient  to  use  several  such  "local"  frames  with  different 
orientations.  Such  frames  are  introduced  to  make  derivations 
shorter,  because  the  probing  wave  scattering  in  each  local  frame 
may  be  treated  as  "almost  forward"  small-angle  scattering.  The 
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lower  boundary  of  the  layer  of  irregularities  is  zd,  the  upper  - 
zu  in  the  global  frame  r.  After  performing  the  Fresnel  expansions 
of  exponents  in  the  (5)  for  small-angle  scattering,  we  obtain  the 
integral  equation  for  U  (7).  This  equation  is  written  in  the  local 
reference  frame,  S(Z-ZQ)  =  (Z'-Zq)T?  +  (Z-Z')R  ;  D(Z-Zq)  = 

(Z’-ZQ)  (Z-Z* ) .  For  making  formulas  shorter  we  assume,  that  the 
transmitter  is  located  exactly  at  the  local  frame  origin,  i.e.  Rq 
=  (P0,Z0)=(0,0,0). 

To  reduce  calculations  and  to  use  known  results  we  transform 
(7)  to  the  parabolic  equation  in  new  variables  £  =  1/Z,  o  =  P/Z, 
$•  =  i/z*.  Then  (S-P' )2/ V  =  (o-S'P'  )2/{V~l)  and  U  satisfies  the 
differential  equation 
d 

(_2ik  —  +  An  -  r2q)  U(?,o)  =  0  (24) 

Equation  (24)  is  derived  from  (7)  by  differentiating  and  making 
use  of  the  following  relation  for  the  fundamental  solution  of  the 
Schrodinger  operator  [24] .  Having  derived  the  parabolic  equation 
we  use  the  methods  developed  to  obtain  equations  for  the  first  and 
the  second  moments  of  U  [26,273.  Let  us  assume  also  that  the 
random  field  q  to  be  Gaussian  and  5-correlated  along  Z, 

BC^  Jiz)=  <  q(P1  ,Z1  )  q(P2,Z2)>  =  6 (Z^Z^  )  A q(?2-P1  )  (25) 

In  the  similar  way  in  the  variables  (a,?;)  it  is  possible  to  derive 
from  (24)  equations  for  the  second  moments  of  the  normalized  field 
7(7?)  =  E(7?)/<E(7?)>, 

F1  Q(P,Z)  =  <t(P2,Z)  qtP1,Z)>  -  the  first  coherence  function, 

-4  *4 

r9  Q(P,Z)  =  <7(P2,Z)  7(P1,Z)>  -  the  second  coherence  function  of 
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the  second  order.  As  the  calculations  are  described  in  the  papers 

mentioned  above,  we  present  only  the  final  result  in  the  variables 

P+  =  (£>i  +^2  )/2,  P  =  311(1  lor  statistically  homogeneous  layer, 

when  the  dependence  on  P+  is  absent. 

-  (26) 
a  V  d  i 

[ik(  +  -  — )  +  Ax  +  Aq(^)  1  I^q  ~  ^2 ,0 —  ^ 

dZ  Z  dP  4k 


The  integral  equation,  corresponding  to  (21 ) ,  has  the  form 


->  i  A  (£’)  ik 

T0  n{P,Z)  =  1  +  —  J  d37?’  — —  T2'Q(P',Z')  exp (W($-P')2) 


'2»°  '  I6iuk  J  V 

The  following  formula  for  T1  1  is  well  known  [26,27]. 


(27) 


1  z  Z' 

Tj^(P,Z)  =  exp  (^2  J  Aq<P“F~>  )  (28) 

Zu 

Eguation  (27)  follows  from  (26).  The  resulting  integral  equation 
is  valid  also  at  PQ  ^  0,  ZQ  ^  0. 

The  integral  equation  (27)  is  the  basis  of  the  statistical 
RT.  According  to  (20)  A  (P)  is  the  correlation  function 
projection,  J*  B(P,Z)  d Z  =  k^(P).  Measured  wave  coherence  functions 
allow  us  to  determine  projections  of  the  complex  potential  q 
correlation  function  by  the  (27,28).  To  begin  with,  let’s  assume, 
that  irregularities  occupy  a  sufficiently  extensive  layer,  oblique 
to  the  probing  wave  propagation  direction.  To  reconstruct 
projection  ^(P)  it’s  necessary  to  determine  the  layer  coordinate 
Zs.  It  may  be  done  by  a  single  receiver.  The  special  procedure  for 
determination  of  the  scattering  layer  Z-coordinate  by  equation 
(22)  was  developed  [1,243. 

In  many  situations  randomly  inhomogeneous  media  represent 


statistically  quasi-homogeneous  fields  with  slowly  varying 
statistical  properties,  i.e.  such  fields,  for  which  the 
correlation  radius  of  the  differential  argument  is  essentially 
less  than  the  scale  of  the  variance  a.  In  other  words,  the 
correlation  function  is  represented  as  a  product  B  (ft^  »ft2 )  = 
o2(ft)/C(Aft,ft),  where  ft  =  (ft.,+  ft2)/2,  Aft  -  ft^-ft^.  A  randomly 
inhomogeneous  field  with  the  costant  correlation  coefficient  £(Ar) 
but  with  varying  fluctuations  variance  a2  (ft)  may  be  called  an 
"additive"  field,  because  the  electron  density  changings  influence 
only  the  fluctuations  intensity,  but  the  correlation  coefficient 

•4 

dependence  on  the  summarize  argument  ft  is  absent.  Assuming  the 
field  q  to  be  8-correlated,  as  was  done  above,  it  is  possible  to 
obtain  similar  equations  by  the  same  technique.  The  correlation 
function  projection  Aq(P)  in  the  integrands  should  be  replaced  by 
o2(Z)  K(P).  Here  K(P)  =  X  K(Aft)  d(AZ’ )  is  the  correlation  function 
projection.  Then  formula  for  r1>t  and  the  equation  for  r2>0  have 
the  form 


r2,0<^>=  1 


+  f  a3»' 

16iuk  J 


o2^’  )K(£') 
V 


(29) 


1  ^  Z' 

r  1#1(^.2T)  =  exp(  — ^  J  o2(Z')  K(p  2~)  d Z')  (30) 

zu 

The  equations  (29,30)  contain  the  unknown  function  a2 (ft).  So,  the 
problem  of  the  correlation  function  reconstruction  is  divided  into 
two  ones:  the  problem  of  reconstruction  of  the  fluctuations 
variance  a2 (ft),  and  the  problem  of  the  K(tti)  reconstruction  given 
the  measured  projections  K(P). 

The  formula  (30)  allows  us  to  determine  integrals  of  o2(Z) 
along  the  direction  ( Z  axis)  of  the  probing  wave  propagation.  If  a 


set  of  different  integrals  is  given  for  a  region  filled  with 
irregularities  we  arrive  at  a  tomography  problem.  Usually  in  the 
experiments  there  exist  a  natural  limit  to  the  number  of 
projections  and  angular  range  of  probing  waves  propagation 
directions.  Hence,  the  problem  when  data  are  known  in  the  small 
angle  range  should  be  solved.  A  similar  problem  was  solved  in  the 
ray  RT.  Information  about  the  fluctuations  intensity  a2  (ft) 
distribution  is  sufficient  to  reconstruct  a  set  of  correlation 
function  projections  given  by  measured  coherence  functions. 
Methods  of  the  statistical  RT  problem  solution  are  similar  to  the 
discussed  above,  with  the  substitution  of  the  product  (^(ftJKtft) 
of  the  known  function  a2  (ft)  and  the  correlation  function 
projection  K(ft)  instead  of  Aq.  A  set  of  projections  K(ft)  allows  us 
to  perform  tomography  reconstruction  of  the  correlation 
coefficient  £(Aft)  or  its  spectrum. 

With  the  help  of  examples  we  want  to  illustrate  the  results 
of  programs  and  numerical  pattering.  Reconstruction  of  electron 
concentration  fluctuations  spectra  was  made  in  two  steps: 
pattering  of  environment  fluctuation  spectra  and  recovery  of  the 
power  angle  spectra  of  the  dispersive  environment. 

I.  Pattering  of  environment  fluctuation  spectra.  In  the  given 
program  <ST-DIR.F0R>  the  account  of  accidental  realisations  of 
electron  concentration  power  fluctuation  spectra  is  made.  The 
pattering  scheme  is  given  of  tabl.3. 

1.  The  sum  of  three  (of  five)  Gauss  lines  with  defferent  width 
and  intensivity  was  taken  as  a  pattern  of  spectra  with  a  few 
maximums.  The  location  of  local  maximums  and  their  width  depend  on 
the  dimentions  of  the  area  (32*32,64*64, etc. )  and  were  given  by 
the  quantity  of  accounts  from  the  beginning  of  the  area. The 
intensivity  of  local  maximums  was  given  in  the  percent  relativity 


from  the  intensivity  of  the  center  maximum. 

2.  The  computation  of  correlation  function  projections  of 
dispersed  field  is  made  in  polar  coordinates  for  the  given 
reception  angles  in  condition  with  equation  (28).  The  projections 
of  environment  correlation  functions  R(p,cJ>n)  calculated  with  the 
help  of  fast  Fourier  transformation  from  the  spectra  (DCp.cJ^)  for 
every  given  reception  angle  in  area  iOs%/2]. 

3.  The  form  of  casual  realisations  of  the  dispersed  field  is 
made  in  two  steps.  At  first,  with  the  help  of  linear  congruent 
method  [28,29]  statisticaly  independent  one-dimentional 
realisations  of  casual  numbers  with  proportional  spectra  were 
obtained.  Then,  with  the  help  of  spectra  method  of  reorganization 
[30]  the  prosess  with  given  correlation  function  ,r1 1  (p.cp^)  was 
formed. 

4.  The  selected  (casual)  correlation  functions  of  the  field 

fSf 

r11  (P»4|q).  were  calculated  from  alone  realisations  of  the  field. 

5.  The  selected  projections  of  correlation  functions  of  the 

fV 

environment  R(p,<f>n)  were  calculated  by  logariphmation  of  the 
function  r11(p,<|)n).  The  appropriate  R ( p , <pn )  casual  power  spectra 

ro 

of  environment  fluctuations  ©(p.c^)  are  defined  with  the  help  of 
standard  procedure  fast  Fourier  transformation.  The  given  model  of 
environment  spectra  (in  Decart  coordinates)  is  given  as  the 
illustration  on  fig.20,a.  Area  dimentions  are  64*64,  width  of 
central  and  local  maximums  on  the  level  of  0.5  their  intensivity 
is  1/8  from  linear  dimentions  of  the  area.  The  local  maximum 
intensivety  is  I-|=0.3IC  (IQ  means  intensivity  of  central  maximum). 
From  the  fig  20, b  one  can  see  casual  spectra  of  environment 
fluctuations,  formed  from  one  field  realisation  for  every  given 
reception  angle.  It  is  easy  to  see,  that  for  more  calculations 
reconstraction  of  angle  spectra  )  the  same  casual  function  need 


smoothing,  l.e.  period  of  Its  oscilations  is  equal  to  the  account 
step  [293.  On  fig.21  statistically  middled  (by  50(a)  ond  100(b) 
realisations)  angle  spectra  of  the  intensivity  of  environment 
fluctuations  ,  took  by  the  given  below  scheme,  are  given. 

II.  Recovery  of  the  power  angle  spectra  of  the  dispersive 
environment . 

In  program  <ST-INV.F0R>  the  account  of  reconstructions  of 
environment  power  angle  spectra  from  the  information  in  finit 
reception  angle  area.  For  solving  the  problem  of  recovery  the  j 

method  of  function  decomposition  into  Lejandr  polynoms  is  used. 

The  system  of  linear  equations,  made  by  this  method,  is  solving  by 
Gordon  method  [303.  On  the  last  step  of  the  account  the  linear 
interpolation  of  the  function  $(p,4>n)  in  Decart  coordinates  is 
made. 

On  fig.  23  one  can  see  reconstructions  of  power  spectra  of 
plasma  fluctuations  calculated  from  one  realisation  of  dispersed 
field  in  each  direction  of  the  reciption  of  the  signal.  The 
account  was  made  for  three-racurs  reception  system  (fig.22,b).  It 
was  supposed  that  angle  areas  makes  0-20°, 35-55°, 65-85°  from 
vertical.  On  fig.23  the  model  and  reconstruction  of  its  casual 
realisation  (the  width  of  central  and  local  maximums  is  1/8  and 
1/10  from  the  linear  dimentions  of  area  (64*64), intensivity 
In=0.3I„)  is  shown.  On  fig. 24  the  results  of  the  same  account  for 
the  dimentions  of  local  and  central  maximums  of  1/16  and  1/8  is 
shown.  Analising  the  results,  we  can  see  that  the  best 
reconstructions  are  obtained  for  the  models  of  sufficiently  narrow 
spectra  (fig.25-28).  In  this  case  (fig.28)  the  recovery  of 
realisation  of  pattern  from  5  local  maximums  (width  is  1/20, 
relativity  of  intensivity  is  1:0. 3:0. 2)  is  sufficiently  good.  On 
fig. 29  the  reconstruction  of  realisation  of  pattern  from  5  local 


maximums  (width  1/16,  relativity  of  intensivity  is  1:0.3:0.2>  is 
given. The  results  of  the  accounts  also  show  that  weak  intensive 
local  maximums  (1:0.1)  recover  sufficiently  good  (fig.  30).  But 
the  essential  part  is  played  by  the  width  of  the  spectra. On  fig. 31 
the  reconstruction  of  model,  which  appropriate  to  the  experimental 
observations [24].  It  is  clear,  that  recovery  of  the  spectra  from 
one  realisation  allow  to  underline  the  particularities  of  its 
structure.  Comparison  analisys  shows  that  in  spectra 
reconstructions  the  locations  of  local  maximums  coincide  with  the 
original  but  the  intensivities  are  differ.  It  future  thes 
difference  can  be  overcomed  by  using  of  correlative  filters. 

The  quality  of  reconstruction  very  much  depend  on  the 
quantity  and  width  of  reception  areas.  On  fig.32,a  reconstruction 
of  model  from  fig.24,  made  by  two-racurs  scheme (fig. 22),  is  given. 
Reception  angles  area  -  0^30°  and  50^80°. ’The  same  reconstructions 
are  made  for  the  experimental  model  (fig. 30).  It  is  clear,  that 
3-racurs  sceme  is  prefered.  On  fig. 33  are  given  spectra 
reconstructions  for  three  reception  areas:  0Q10°, 40^50°, 1 70^80°, 
which  are  in  connection  with  fig. 32.  It  is  clear  that  to  give  more 
quality  to  the  reconstruction  the  width  of  reception  area  should 

have  optimum  value. 


The  description  of  Fortran  Programs 


System  Requirements 


-  Computer: 

-  Operating  System: 

-  Memory: 


-  Hard  Disk  Space: 

-  Software: 


IBM  AT  or  compatible  (with  coprocessor) 
MS-DOS  or  PC-DOS  version  3.0  and  later 
at  least  Extended  memory  4  Mbytes 
(depends  on  size  of  frame  and 
namber  of  rays) 

4  Mbytes 

compiler  1.4e  and  linker  2. 2d 

NDP -FORTRAN-386 (c )  MicroWay  or  later 


1.  Program  <ST-DIR.for> 

This  program  solves  direct  problem,  namely, 
determines  the  model  structure  and  calculates 
the  random  realizations  of  angle  power  spectrum. 


Input  parameters  and  files: 

M  -  size  of  frame 
NM  -  size  of  half  frame 
JMN  -  number  of  rays  in  [0;tu/2]  range 
KR1  -  local  maximum  size 
KR  -  central  maximum  size 
N00  -  spectrum  coordinate  of  local  maximum 
OTT  -  level  of  local  maximum  power 
array  HH  -  model  of  angle  power  spectrum  in  polar 
coordinates 

array  SPAD  -  correlation  function  model  of  scattering  field 
Subroutine  SMAXMI  calculates  maximum  and  minimum  of 

one-dimension  function 

Subroutine  SMAXNI  calculates  maximum  and  minimum  of 

two-dimension  function 

Subroutine  SPGEN  calculates  balance  multipliers 
Subroutine  RANFII  calculates  random  fields  (complex) 
Subroutine  CORPAT  calculates  random  correlation  function 
Subroutine  F0RT1  calculates  one-dimension  Fast  Fourier 

Transformation 

Function  RANI  calculates  uniform  random  number 


Output  files: 

<RANDSP.GRD>  -  contains  random  realization  of  angle 
power  spectrum,  size  (JMN:M) 

Execution 

177  st-dir.for 

RUN 

st-dir.exp 

3.  Program  <ST-INV.for> 

This  program  solves  inverse  problem  for  statistical  RT, 
namely,  calculates  reconstructions  of  power  spectrum  of  electron 
density  in  Decart  coordinates  for  varios  namber  of 
receivers. 

Input  parameters  and  files: 

Parameters 
M  -  size  of  frame 
NM  -  half  size  of  frame 
NN  -  half  size  of  frame  +1 
ITER  -  number  of  iterations 
JPI  -  number  of  reseivers 
JMN  -  numer  of  rays  in  LQ,%/2]  range 
JT1 ,  JT2,  JT3,  JT4, JT5, JT6  -  boundaries  of  angle  ranges 

of  reception  in  degrees 

NU  -  number  of  rays  in  angle  ranges  of  reception 
NXX  -number  of  lejandr  cofficient 

PILE  <RANDSP.GRD>  contains  realization  of  angle  power  spectrum 

(OUTPUT  file  of  progpam  ST-DIR.for) 

OUTPUT  FILE: 

<  RECONS. GRD>  contains  the  reconstraction  of  realizations 
angle  power  spectrum  in  Decart 
coordinates , s ize (M: M ) . 

Execution 

f77  st-inv.for 

RUN 

st-inv.exp 


% 


References. 


1.  V.E. Kunitsyn,  E.D. Tereshchenko,  Tomography  of  the  Ionosphere, 

Moscow,  Nauka,  1991  (in  Russian). 

2.  V.E. Kunitsyn,  "Inverse  Scattering  Problem  in  the  Layered  Medium", 

Moscow  University  Physics  Bulletin,  Vol.26,  No  6,  pp. 33-40,  1985. 

3.  V.E. Kunitsyn,  "Determination  of  the  structure  of  ionospheric 

irregularities".  Geomagnetism  and  Aeronomy,  Vol.26,  pp. 75-81 ,  1986 

4.  V.E. Kunitsyn,  "Diffraction  Tomography  in  the  Fresnel  Zone", 

Moscow  University  Physics  Bulletin,  Vol.41 ,  No  2,  pp. 43-48,  1986. 

5.  V.E. Kunitsyn,  N.G. Preobrazhensky,  E.D. Tereshchenko, 

"Reconstruction  of  ionospheric  irregularities  structure  based 
on  the  radioprobing  data",  Dokl.Akad.Nauk  USSR,  Vol.306, 
pp. 575-579,  1989. 

6.  V.E. Kunitsyn,  E.D. Tereshchenko,  The  reconstruction  of  the 

ionosphere  irregularities  structure ,  Preprint  Polar  Geophys. 

Inst.,  90-01-69,  pp.1-60,  1990. 

7.  V.E. Kunitsyn,  "Diffraction  Tomography  based  on  small-angle 

scattering  data".  Transactions  SPIE,  Vol.1843,  pp. 172-1 82,  1992. 

8.  V.E. Kunitsyn,  E.D. Tereshchenko,  "Radio  Tomography  of  the 

Ionosphere",  Antennas  Rropagation  Magazine,  Vol.34,  pp. 22-32,  1992 

9.  V.L. Ginzburg,  Propagation  of  electromagnetic  waves  in 

plasma.  New  York,  Gordon  and  Breach  Science  Publishers,  Inc.,  1961 

10.  L.D. Landau,  E.M.Lifshits,  Statistical  physics ,  New  York, 

Pergamon  Press,  1980. 

11.  L.D.Faddeev,  "Three-dimensional  inverse  scattering  problem", 

J.Sov.Math. ,  Vol  5,  pp.  334-360,  1976. 

12.  R.G. Newton,  "Inverse  scattering  II.  Three  dimentions"  J.Math. 

Phys.,  Vol. 21,  pp. 1698-1715,  1980. 

13.  K.Chadan,  P.G. Sabatier,  Inverse  problems  in  quantum  theory , 

New  York,  Springer-Verlag,  1977. 

14.  A. J.Devaney,  "Geophys.  diffraction  tomography",  IEEE  Trans. 

Geosci.  and  Remote  Sencing,  Vol. 22,  pp.3-13,  1984. 

15.  R.K. Mueller,  M.Kaveh,  G.Wade,  "Reconstructive  tomography  and 

applications  to  ultrasonics",  Proc.  IEEE,  Vol. 67,  pp. 567-587,  1979 

16.  D.Lesselier,  D.Vuillet-Laurent,  F.Jouvie  at  al.,  "Iterative 

solutions  of  some  direct  and  inverse  problems  in  electromagnetics 
and  acoustics".  Electromagnetics,  Vol. 5,  pp. 147-189,  1985. 


17.  K.  J.Langeberg,  "Applied,  inverse  problems  for  acoustic, 

electromagnetic  and  elastic  wave  scattering",  in  P. Sabatier  (ed. ), 
Basic  methods  of  tomography  and  inverse  problems,  Bristol, 

Adam  Hilger,  pp. 127-467,  1987. 

18.  W.Tabbara,  B.Duchene,  G.Pichot.  et  al,  "Diffraction  tomography: 

contribution  to  the  analysis  of  some  applications  in  microwaves 
and  ultrasonics".  Inverse  Problems,  Vol.4,  pp. 305-331 ,1988. 

19.  W.K. Pratt,  Digital  Image  Processing,  New  York,  A  Wiley 

Interscience  Publ.,  1978. 

20.  D.E.Dadgeon,  R.M.Mersereau,  Multidimensional  Digital  Signal 

Processing,  Prentice-Hall,  Inc.,  1984. 

21.  D.Slepian,  "About  bandwidth",  Proc.IEEE,  Vol.64,  No. 3,  pp.292-300. 

22.  H.J. Landau,  H.O. Poliak,  "Prolate  spheroidal  wave  functions, 

Fourier  analysis  and  uncertainty".  Bell  Syst.Tech. J. ,  Vol.40, 
pp.  65-84;  Vol.41.  pp. 1295-1 336,  1961. 

23.  E.S. Andreeva,  V.E. Kunitsyn,  "The  reconstruction  of  the  projections 

of  scattered  irregularities  by  inesact  wavefield  data".  In: 

Waves  propagation  and  diffraction  in  inhomogeneous  media, 

Moscow,  Moscow  Inst.  Physics  and  Technology,  pp.  9-14,  1989. 

24.  V.E. Kunitsyn,  E.D. Tereshchenko,  "Determination  of  the 

turbulent  spectrum  In  the  ionosphere  by  a  tomographic  method", 
Joum. Atm. Terr. Phy s. ,  Vol.54,  pp. 1275-1282,  1992. 

25.  A.V.Galinov,  V.E. Kunitsyn,  E.D. Tereshchenko,  "A  tomographic 

approach  to  the  investigation  of  a  randomly  nonuniform 
ionosphere".  Geomagnetism  and  Aeronomy,  Vol.31,  No. 3, 
pp. 351 -356,  1991 

26.  V.I.Tatarski,  Wave  propagation  in  a  turbulent  medium. 

New  York,  McGraw-Hill,  1961. 

27.  A.Ishimaru,  Wave  propagation  and  scattering  in  random  madia, 

Orlando  (Florida),  Academic  Press i  1978. 

28.  I.M. Sobol,  Monte-Carlo  method ,  Moscow,  Nauka,  1978  (in  Russian). 

29.  S.Chakrabarti,  S.K.Mitra,  "Calculations  of  two  dimensional 

digital  filters  by  spectral  convolution",  Proc.IEEE, 

Vol.65,  No. 6,  1977. 

30.  R. Gordon,  R. Bender,  A. Lent,  Joum. Theor. Biol. ,  Vol.29,  1970. 


A  (KM) 


Pc 

pi2 

0.154 

0.183 

0.345 

0.436 

0.439 

0.675 

0.672 

0.882 

0.093 

0.101 

0.241 

0.253 

0.468 

0.475 

0.654 

0.717 

1 


Pattern  of  angle  spectrum 

®<p.V 


TaM.  3 


•ST 


Ff2, 


Development  of  Algorithms  and  Programs  for  the 
Diffraction  and  Statistical  Radiotomography  of  the  Ionosphere. 

Final  (12-month.)  report  (EOARD  Contract  SPC  93-4063) 


CONTENT  PP- 

1 .  Introduction.  1 

2.  Physical  and  mathematical  formulation  of  the  radio  tomographic 

sounding  problem.  4 

3.  The  solution  of  the  diffraction  radiotomography  problem 

(inverse  scattering  problem).  9 

Description  of  the  program  for  calculation  of  scattered 


radiowaves  fields  (DT  DIR. FOR).  30 

Description  of  the  program  for  solution  of  the  diffraction 

radiotomography  problem  (DT  INV.FOR).  32 

Description  of  the  program  for  estimation  error  of  RT 

reconstruction  (DT  ERR. FOR).  32 

4.  The  solution  of  the  statistical  radiotomography  problem 

(statistical  inverse  scattering  problem)  34 

Description  of  the  program  for  calculations  of  the  statis¬ 
tical  characteristics  of  scattered  waves  (ST  DIR. FOR).  42 

Description  of  the  program  for  solution  of  the  statistical 

radiotomography  problem  (ST  INV.FOR).  43 

Referencees  44 

Figures  and  tables  46-81 


Principal  investigator 


Professor  V. Kunitsyn 


1 .  Introduction 


The  aim  of  this  work  was  to  develop  efficient  algorithms  and 
programs  for  the  solutions  of  the  ionospheric  diffraction  and 
statistical  radiotomography  (RT)  problems.  The  diffraction  and 
statistical  RT  problems  are  related  with  reconstruction  of  the 
scattering  ionospheric  irregularities  structure.  The  use  of 
tomographic  methods  is  perhaps  an  unavoidable  stage  in  evolution 
of  nearly  all  diagnostic  systems.  At  a  sufficiently  high  level  of 
development  of  remote  sensing  technology  and  data  processing 
resources,  it  become  possible  to  reconstruct  the  spatial  structure 
of  a  medium  on  the  basis  of  tomography.  Tomographic  approaches 
have  already  transformed  sensing  methods  in  many  fields  and 
provided  fundamentally  new  results.  The  major  achievements  of 
tomography  in  medicine  and  molecular  biology  are  widely  known. 
Tomographic  methods  have  made  it  possible  to  detect  previously 
unknown  phenomena  in  geophysics  (seismo tomography  of  the  Earth  and 
acoustic  tomography  of  the  ocean).  At  the  present  stage,  radio 
sounding  technology  makes  it  possible  to  use  satellite  resources 
to  perform  ionospheric  sounding  in  a  wide  range  of  different 
positions  of  the  transmitting  and  receiving  systems  and  to  use 
tomographic  methods.  In  connection  with  this,  work  on  RT  of  the 
ionosphere  has  began  actively  in  recent  years. 

The  purpose  of  this  paper  is  to  describe  shortly  the  RT 
methods  Cl -8]  and  program  developed  for  reconstruction  ionospheric 
irregularities  structures  together  with  some  results  of  computer 
simulation.  We  are  dealing  with  reconstruction  of  the 
irregularities  of  the  electron  density  and  of  the  effective 
collision  frequency  by  scattered  radio  waves  field.  Mathematicians 
used  very  often  designations  the  Inverse  problem  (IP)  or  the 


inverse  scattering  problem  (ISP)  for  such  problems  of  structure 
reconstruction  owing  to  the  scattered  field.  Here  we  used  the 
terms  IP  or  ISP  and  diffraction  tomography  as  synonims.  The 
ionosphere  has  a  rather  complex  structure;  local  irregulariries  of 
various  scales,  including  turbulent  regions,  are  present  in 
addition  to  a  quasi-stratified  background  with  large  scales. 
Therefore,  problems  of  satellite  RT  of  the  ionosphere  should  be 
divided  into  statistical  RT  and  deterministic  RT  problems.  The 
latter,  in  turn,  are  divided  into  problems  of  diffraction  RT  and 
ray  RT  of  large  structures,  where  diffraction  effects  are  not 
significant. 

There  is  no  special  reason  to  provide  a  strict  classification 
of  any  field  of  science  developing,  but  certain  delimitations  and 
explanations  of  the  terminology  must  be  provided  so  that  we  will 
not  later  be  reproached  for  "inventing"  new  terms  unnecessarily. 
In  the  cases  where  any  projections  or  "cross  sections"  of  an 
inhomogeneous  object  (or  a  certain  transformation  of  the  object, 
such  as  Fourier  transform)  are  known  from  remote  sensing  data  and 
the  problem  is  to  reconstruct  the  structure  of  the  object,  the 
problem  should  be  considered  as  a  tomographic  problem.  At  present, 
the  mathematical  foundation  for  tomography  is  related  to  integral 
geometry,  where  it  is  necessary  to  reconstruct  an  object  using 
data  on  It  in  the  form  of  integrals  with  respect  to 
small-dimension  manifolds.  Therefore,  the  term  tomography  is 
understood  not  in  narrow  initial  sense,  as  layer-by-layer  study  of 
the  structure  of  the  inhomogeneous  objects,  but  rather  in  the 
wider  sense,  as  Recording  of  projections  or  cross-sections  of  an 
object  and  subsequent  reconstruction  of  the  structure  of  the 
object  from  them.  The  projections  of  an  object  are  various 
integrals  with  respect  to  small-dimension  manifolds.  The  range  of 


problems  of  radio  remote  sensing  of  the  ionosphere  using 
satellites  examined  in  this  paper  includes  reconstruction  of  the 
structure  of  ionospheric  irregularities  from  files  of  various 
tomographic- type  data  (projections,  cross  sections),  which 
provides  the  basis  for  using  the  term  RT  of  the  ionosphere. 

Problems  of  satellite  RT  should  be  divided  into  deterministic 
and  statistical  problems.  In  the  case  of  a  deterministic  RT 
problem,  it  is  necessary  to  reconstruct  the  structure  of  certain 
large  irregularities  or  a  group  of  irregularities.  If  a  large 
number  of  irregularities  occupies  a  certain  region  in  space,  it  is 
unadvised  to  reconstruct  the  structure  of  all  realizations 
(generations)  of  irregularities  every  moment  of  time.  Here,  it 
makes  sense  to  pose  the  problem  of  reconstruction  of  the 
characteristics  of  full  ensemble  of  irregularities,  the  structure 
of  the  statistical  characteristics  of  the  irregularities  such  as 
the  correlation  function  of  the  electron  density,  etc. 

We  will  single  out  a  number  of  fundamental  features  of 
problems  of  RT  of  the  ionosphere  or  the  IP  of  reconstruction  of 
the  structure  of  the  inhomogeneous  ionosphere.  The  dimensions  of 
the  transmitting  and  receiving  systems  feasible  in  practice  are 
much  less  than  the  distances  from  them  to  the  irregularities  to  be 
reconstructed,  which  are  hundreds  of  kilometers,  i.e.,  the 
aperture  angles  are  small.  Since  it  is  extremely  difficult  and 
expensive  to  create  transmitting  and  receiving  systems  with  a 
large  number  of  receivers,  the  necessity  for  aperture  synthesis  in 
one  coordinate  becomes  clear.  Synthesis  apertures  can  be  realised, 
for  example,  using  a  moving  transmitter  on  a  satellite.  We  will 
emphase  that  satellite  radio  sounding  makes  it  possible  in 
practice  to  obtain  files  of  various  tomographic  data  and  to 
actually  realize  RT  of  the  inhomogeneous  ionosphere.  In  view  of 
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the  small  aperture  angles,  it  is  advisable  to  pose  the  problem  of 
reconstruction  of  the  structure  of  irregularities  with  dimensions 
singficantly  exceeding  the  wavelength.  Therefore,  the  IP  of 
reconstruction  of  the  structure  of  scatterers  which  are  large  in 
compasion  with  the  wavelength  will  be  examined  here,  but  the 
dimensions  of  the  irregularities  may  be  both  greater  than  or  less 
than  Fresnel  zone,  and  diffraction  effects  must  be  considered  in  a 
number  of  cases. 


2.  Physical  and  mathematical  formulation  of  the  radio  tomographic 
sounding  problem. 

Radio  wave  propagation  in  near-earth  space  and  scattering  in 
plasma  irregularities  in  the  ionosphere  are  discribed  by  a  system 
of  Maxwell  equations  with  the  corresponding  material  equations 
[93.  In  the  general  case,  the  dielectric  permittivity  of  the 
plasma  is  a  tensor  value  due  to  spatial  dispersion,  even  in  the 
absence  of  a  magnetic  field.  However,  the  spatial  dispersion  can 
be  disregarded  for  ionospheric  irregularities  of  both  natural  and 
artificial  origin,  since  the  thermal  velocities  of  the  electrons 
in  the  plasma  are  significantly  less  then  than  the  speed  of  light, 
i.e.,  the  approximation  of  a  "cold"  plasma  is  valid.  Likewise,  the 
velocities  of  diffusion,  mixing  and  other  transport  processes,  as 
well  as  the  velocity  of  the  transmitting  and  receiving  systems, 
does  not  exceed  10  km/s  (v/c  <  3*10-4),  therefore  ,  the 
quasistationary  approximation  is  valid;  within  the  framework  of 
it,  the  "slow"  time  dependence  of  the  characteristics  of  the 
medium  and  fields  can  be  considered  as  a  dependence  on  the 
parameter. 


In  place  of  the  equations  for  the  field  £(r,t)  of  a  radio 
wave,  it  is  convenient  to  use  equations  for  the  complex  amplitudes 
1  of  the  corresponding  monochromatic  components,  i.e.,  it  is 
advisable  to  introduce 

2(r,t)  =  S(r,t)e_iut  +  3*(r,t)eiwt. 

For  brevity,  the  complex  amplitudes  S(r,t)  with  a  "slow"  time  will 
be  called  the  field  hereafter.  Subject  to  the  observations  made 
above,  within  the  framework  of  the  quasistationary  approximation 
and  a  "cold"  plasma,  from  the  Maxwell  equations  we  have  the 
equation  for  the  field 


AE  + 


i  E  -  grad  div  E  =  0, 


where  s  is  the  complex  dielectric  permittivity  tensor,  f  =w/: 2%  is 
the  frequency.  Analysis  of  radio  wave  propagation  in  an 

A. 

inhomogeneous  magnetoactive  plasma  with  a  tensor  dependence  of  s 

is  a  difficult  problem.  There  are  no  algorithms  for  calculating 

the  parameters  of  radio  signals  and  the  characteristics  of  the 

scattered  fields  in  the  general  case  of  an  arbitrary  dependense  of 

s(r)  on  the  coordinates.  Therefore,  the  statement  of  the  ISP  of 

reconstruction  of  the  structure  of  the  tensor  characterising  an 

ionospheric  irregularity  using  data  on  the  measured  field  is 

unrealistic  at  present,  and  hardly  advisable  in  general.  The 

non-diagonal  elements  of  the  tensor  can  be  disregarded  at  high 

frequencies,  since  they  do  not  exceed  the  square  of  the  ratio  of 

the  plasma  frequency  fN  to  the  sounding  frequency,  i.e.,  ~ 

(fjj/f )2(fH/f )  (fH  is  the  gyrofrequency)  [93.  For  example,  the 

6  —3 

typical  maximum  concentration  in  the  ionosphere  is  ™  10  cm  . 
For  it  and  where  f  >  50  Mhs,  (f^/f  )2(fH/f  >  ^  0,001.  Likewise,  the 
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last,  "depolarisation"  term  of  equation  (1),  the  order  of  which  is 
determined  by  the  ratio  of  the  emission  wavelength  to  the 

characteristic  scale  of  variation  of  the  concentration,  can  also 
be  disregarded  at  high  frequencies. 

Thus,  in  the  case  of  high  sounding  frequencies  vector 

equation  (1 )  decomposes  into  three  scalar  equations,  and  it  is 
sufficient  to  examine  the  equation  for  one  field  component 

A  E  +  k2s(r,k)  E  =  0,  (2) 

4,nrJJ(r) 

where  a(r,u))  =  1 - 5 - — - »  k  =  2Tuf/c  =  w/c  is  the 

k2(1  +  iv(r)/o)) 

wave  number,  r0  is  the  classical  electron  radius,  v(r)  is 
effective  collision  frequency,  N(r)  is  the  electron  concentration; 
"ion"  contribution  to  the  dielectric  permittivity  can  be 
disregarded  [33.  The  relation  for  e  has  the  same  form  in  both  the 
International  System  of  Units  (SI)  and  CGS  system;  only  the 
expression  for  rg  changes. 

Then,  it  will  be  convenient  to  introduce  the  complex  function 
q(r,to)  =  k2(1-e)=  4TCreN(r)/(1  +  iv(r)/w)  and  divide  it  into  two 
parts:  q  -*  qQ(z,u))  +  q(r,w),  one  of  which  corresponds  to  the 
regular  stratified  ionosphere  with  the  dependencies  NQ(z)  and 
v  (z),  and  the  other  of  which  is  related  to  the  three-dimension 
(3D)  irregularities  on  the  background  of  the  stratified  medium 
N(r)  and  v(r).  The  specific  and  concrete  form  of  relations  for 
vQ(z)  and  v(r)  depends  on  the  specifics  of  the  problem  [13.  By 
reconstructing  the  structure  of  the  complex  function  q(r,w),  it  is 
also  possible  to  reconstruct  the  structure  of  these  3D 
disturbances  with  known  information  on  the  regular  large  scale 
ionosphere.  Therefore  we  will  deal  with  the  two  problems:  first, 
the  problem  of  reconstruction  of  the  regular  large  scale 


ionosphere  q0(z,w)  -  ray  RT  problem;  second,  the  problem  of 
reconstruction  the  2D  or  3D  irregularities  structure  - 
dif fractional  and  (or)  statistical  RT  problem.  The  introduction  of 
q  is  also  justified  by  the  fact  that  the  quantity  q  ~  (1  -  e) 
applies  to  the  generalized  susceptibilities  [10]  characterizing 
the  response  of  the  medium  to  the  field.  Furthermore,  the  scalar 
Helmholtz  eguation  assumes  the  form  of  stationary  Schrodinger 
equation,  where  q  is  the  complex  potential: 

AE  +  k2E  -  q0(z,w)E  -  q(r,o))E  =  S(r  -  rQ)  (3) 

The  delta  function  in  the  right  side  of  (3)  corresponds  to  the 
point  sounding  source  in  the  majority  of  satellite  radio  sounding 
experiments,  since  the  wave  can  be  assumed  to  be  spherical  within 
the  major  lobe  of  the  directional  pattern  of  an  arbitrary  source. 
The  dimensional  constant  AQ  in  front  of  the  5-function  in  equation 
(3)  is  omitted  so  that  it  need  not  be  rewritten  in  all  subsequent 
formulas.  The  constant  is  expressed  by  means  of  the  total  power  of 
a  point  source  [1j.  Here,  A0  corresponds  to  the  amplitude  value  of 
the  field  of  a  point  source,  i.e., 

EQ(r)  =  -  Aq  exp(ikr)/(4icr). 

Hereafter,  we  can  always  return  to  the  dimensional  field  E  -*  AQE 
when  necessary. 

The  ISP  for  equation  (3)  can  be  formulated  in  the  following 
manner:  it  is  necessary  to  reconstruct  the  structure  of  the 
irregularities  from  measurements  of  the  field  in  a  certain  limited 
region  of  a  fixed  surface  in  a  limited  range  of  frequencies  and 
positions  of  the  sounding  source.  Knowledge  of  the  complex 
function  q(r,(jj)  makes  it  possible  to  reconstruct  the  structure  of 
both  N(r)  and  v(r).  Here  q(r,w)  is  a  finite  function,  since  the 
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irregularities  of  N(r)  are  confined  in  space.  The  regular 
stratified  ionosphere  is  assumed  to  be  known  here,  i.e.,  the  IP  of 
reconstruction  of  the  regular  ionosphere  is  solved  by  other, 
sufficiently  well-developed  methods,  or  its  influence  can  be 
disregarded,  for  example,  as  in  the  UHF/VHF  waveband. 

The  sheme  of  the  experiments  on  diffraction  and  statistical 
RT  is  presented  in  Pig-1.  The  satellite,  with  a  transmitter  on 
board  moves  at  altitude  z  =zQ;  the  receiving  system  is  located  in 
the  plane  z  =  zR  (on  the  ground).  The  receiving  system  can  be  a 
transverse  array  of  receivers  for  diffraction  RT  or  statistical  RT 
(or  a  collection  of  transverse  arrays).  Receivers  of  array 
arranged  along  (the  y-axis)  a  line  transverse  to  the  satellite 
path  (the  x-axis),  arrays  being  separated  by  distances  of  hundreds 
of  kilometers.  The  scatterer  (a'  group  of  irregularities)  is 
located  near  the  altitude  z  =  z  . 

S 


3.  The  solution  of  the  diffraction  radiotomography  problem 
(inverse  scattering  problem). 

Rigorous  results  of  solution  of  IP  (3)  for  the  complex 
potential  do  not  exist.  We  will  emphasize  that  in  the  case  of 
propagation  of  waves  of  various  nature  in  the  actual  physical 
media  the  complex  nature  of  the  potential  q  has  a  fundamental 
character  and  the  corresponding  dispersion  relations  relating  the 
real  and  imaginary  parts  of  the  potential  can  be  derived  [103.  The 
fundamental  results  on  solution  of  the  3D  inverse  problem  for  the 
Schrodinger  equation  with  the  real  potential  q(r)  not  depend  on 
the  frequency  (where  q0(r,w)  =  0)  obtained  in  [11,123  should  be 
noted  here.  As  to  the  complex  nature  of  the  potentials,  there  no 
rigorous  results  for  the  ID  IP,  i.e  the  extension  of  the  familiar 
Gelfand-Levitan-Marchenko  algorithm  to  this  case  is  unlmown,  even 
more  so  for  a  arbitrary  dependence  of  the  potential  from  the 
energy  (frequency  co)  [133. 

Approximate  approaches  to  solution  of  the  3D  ISP  are 
therefore  of  undoubted  interest.  The  approximate  approach  to 
solution  of  inverse  RT  problem  is  justified  by  the  properties  of 
the  regular  ionosphere,  the  typical  scales  of  essential 
disturbances  of  which  significantly  exceed  the  radio  sounding 
wavelength.  Therefore,  the  geometrical  optics  approximation  is 
applicable  to  description  of  wave  propagation  on  the  background  of 
the  regular  ionosphere.  The  scattering  in  weak  irregularities  is 
described  within  the  framework  of  such  familiar  approximations  as 
the  Bom  approximation  and  Rytov  approximation  (BA  and  RA).  More 
specialized  methods  iiave  been  used  to  calculate  the  scattering  in 
strong-scattering  Irregularities. 

There  are  many  papers  connected  with  the  approximate 


dt 


approaches  to  ISP  [14-183.  The  scattering  by  weak  irregularities 
was  described  within  the  BA  and  RA.  The  iterative  methods  were 
used  for  strong-scattering  irregularities.  But  there  are  some 
specific  features  of  the  ionosphere  problems,  namely:  the  size  of 
irregularities  is  very  big  (it  can  exceed  the  probing  wavelength 
in  103-  104  times),  the  ionosphere  irregularities  can  be  either 
weak-scattering  or  strong-scattering.  It  is  very  difficult  to  use 
iterative  methods  for  big  size  strong-scattering  irregularities 
directly.  However  it  is  sufficient  to  use  the  small-angle 
approximation  according  to  the  big  size  of  irregularities  or  their 
detailes.  In  what  follows  we  consider  the  methods  of  ISP  solutions 
only,  based  on  the  asymptotic  approximation  for  the  forward 
small-angle  scattering.  It  should  be  noted  it  is  particularly 
interesting  for  tomographic  methods  to  consider  such  cases  of 
small-angle  scattering,  when  a  wavelength  of  sounding  radiation  is 
much  less  compared  with  the  scales  of  object  specific  detailes. 
This  is  necessary  condition  to  reconstruct  the  complicated 
internal  structure  of  the  object. 

In  the  beginning  it  is  useful  to  consider  the  ISP  In  the 
high-frequency  limit  i.e.  the  source  field  freguency  is 
essentially  higher  than  the  critical  frequency  of  the  regular 
ionosphere.  In  this  case  the  refraction  index  is  close  to  unity, 
and  the  influence  of  the  regular  Ionosphere  may  be  neglected. 
Then,  Instead  of  (1),  the  following  equation  should  be  used: 

AE  +  k2E  -  q(r,w)E  =  6(r-rQ)  (4) 

This  differential  equation  is  equivalent  to  the  Lippmann-Schwinger 
integral  equation 

E(r)  =  E0(r)  +  J  G(r-r')  E(r* )  q(rf,a))  d3r’ 


(5) 


where  G(r)=-(4icr)  1exp(ikr)  is  the  Green’s  function  for  a  vacuum 
and  EQ(r)=G(r-r0)  is  the  source  field. 

It  is  necessary  to  bear  in  mind  that  the  potential  q  depends 
on  N  and  v  as  follows  (the  usual  condition  for  ionospheric 
radioprobing  v  «  w  is  taken  into  account )  E 1 3 

q  =  to  N(r)/(1  +  iv(r)/w)  «  4nxJI(r)(1  -  iv(r)/w)  (6) 

T.  0  C 


For  the  case  of  small-angle  scattering  the  transmitter,  a 
scattering  object  and  the  receiver  are  situated  along  a  straight 
line  approximately.  Assuming  that  the  direction  of  this  straight 
line  is  closed  to  the  direction  of  z-axis  one  can  obtain  an 
approximate  equation,  describing  the  forward  small-angle 
scattering,  instead  of  (5)  El, 7 3 

U(r)  =  1  +  f  F(r,r’ )  U(r’)  q(r’,w)  dV  =  1  +  FUq  (7) 


Here  U=E/E0.is  the  normalized  field.  Under  the  derivation  (7)  the 
Fresnel's  approximation  was  used  for  the  Green  function  G  either 
as  for  the. sounding  wave  field.  After  some  transformations  we  come 
from  (5)  to  the  formula  (7),  where  the  kernel  F  of  the  equation 
(7)  contains  coordinates  of  the  transmitter  Zq,  the  finite 
scatterer  zg  and  the  receiver  z^ 


«| 

F(r,r’ )  - - exp  (i-^Wp'-  s)2), 

4tc£  ^ 


(8) 


®  -  Pr^Z0  Zs>/(Z0  ZR>  +  P0^Zs  ZR^Z0  ZR^: 


C  -  (z0  "  zs)(zs  ~  zr)/(z0  ZR^* 

In  general  case  the  formulae  for  s  and  C  3iave  to  include  the 
integration  variable  z’  instead  of  zs-  But  if  the  object  scale  is 
less  than  longitudinal  Fresnel  resolution  the  varriable  z'  can  be 


changed  on  the  constant  zs- 

There  is  another  form  of  the  equation  (7),  using  the  complex 
phase  $  =  In  E,  ©  -  ©Q  =  In  U 

©  -  ©0  =  In  (1  +  Fq  exp(©  -  ©0))  (9) 

which  includes  the  RA  (the  first  iteration  of  the  non-linear 
equation  (9)) 

©  -  ©0  =  In  (1  +  Fq)  “  Fq  or  U  =  exp(Fq)  (10) 

Taking  into  account  (7,10)  and  making  use  of  the  Fresnel's 
transformation  we  have  the  algorithm  for  reconstruction  of  2D 
cross-section  of  potential  q  [1,23 

q2(p,  a)  s  Jq(r,G»  dz  =  (k/2iuC)2  J  V(s)  exp(-ik(s-p)2/2C)d2s. 

(11) 

The  integral  projection  q^(p)  can  be  obtained  after  transformation 
of  field  data.  For  case  of  the  BA  we  have  V=-4'iuC(E-E0)/E0  and  for 
the  case  of  the  RA  V  =  -4tc£(©-©q)  [1,23.  Changing  the  direction 
of  the  sounding  wave  incident  ("turning"  the  z-axis)  one  can 
obtain  a  set  of  2D-projections,  which  makes  it  possible  the 
tomography  reconstruction  of  the  3D  structure. 

Determination  of  the  linear  potentials  q  (p,w)  (11)  from 

complex  potential  q(r,io)  in  the  case  of  forward  scattering  reduces 
the  ISP  to  a  problem  of  tomographic  reconstruction  -  the  problem 
of  reconstruction  an  object  from  projections.  The  intensive 
development  in  the  past  decade  of  tomographic  methods  for  studying 
the  structure  of  objects  is  to  a  significant  extent  the  result  of 
advances  in  x-ray  tomography.  The  reconstruction  algorithms  used 
in  practical  x-ray  tomograpiiy  are  based  on  a  linear  approximation 
of  the  ray  trajectories.  In  mathematical  respects,  such  problems 


reduce  to  reconstruction  of  the  attenuation  function  or  refractive 
index  from  group  of  linear  integrals  and  to  reconstruction  of  an 
object  from  its  small-dimensional  projections.  In  addition  to 
x-rays,  practically  all  known  types  of  emissions  and  waves  are  now 
used  for  tomographic  purposes.  The  linear  approximation  frequently 
does  not  provide  good  results  in  tomographic  investigations  using 
optical  and  ultrasonic  waves,  microwaves  and  other  waves. 
Therefore,  methods  for  reconstruction  with  consideration  of 
reflection  and  diffraction  effects  have  been  developed  intensively 
in  recent  years,  and  a  special  term  -  diffraction  tomography  -  has 
appeared  [ 8-1 2,43. 

The  solutions  of  the  IP  of  reconstruction  of  the  structure  of 
a  scatterer  with  diffraction  considered  in  the  case  of  "forward" 
scattering  obtained  above  can  also  be  applied  to  the  field  of 
diffraction  tomography.  ISP  reduces  to  the  tomographic  problem  of 
reconstructing  a  3D  object  from  2D  projections  qs(p,w)  (11).  By 
rotating  the  object  or  rotating  the  transmitting  and  receiving 
system  in  relation  to  the  object  (axis  z),  it  is  possible  to 
obtain  a  series  of  linear  integrals  of  qz(p,w),  which  can  be  used 
for  tomographic  reconstruction  of  the  refractive  index.  However, 
the  problem  of  reconstructing  the  2D  structure  from  a  set  of 
ID  "cross  sections",  i.e.,  functions  measured  by  one  receiver  as 
the  transmitter  on  a  satellite  moves,  is  also  a  tomographic 
problem.  Furthermore,  the  physical  meaning  of  the  quantities 
measured  by  one  receiver  is  of  integrals  over  the  cross  section  of 
the  spectrum  of  the  irregularity.  A  set  of  such  integrals  makes  it 
possible  to  reconstruct  the  2D  cross  section  of  spectrum,  which  is 
equivalent  to  reconstruction  of  2D  integrated  projection  of 
qi7(p,w).  The  set  of  2D  projections  makes  it  possible  to  perform 
tomographic  reconstruction  of  the  3D  structure. 


The  influence  of  diffraction  reduces  to  the  fact  that  the 
linear  integral  (11)  is  dependent  not  only  on  the  field  "on  the 
ray",  bur  also  on  the  field  in  vicinity  of  the  intersection  point 
of  the  ray  and  recording  plane  (z  =  zR).  If  the  emission  frequency 
is  high,  then  at  the  limit  w  »  the  integral  (11 )  should  only  be 
dependent  on  the  field  on  the  ray  (the  eikonal  approximation  with 
a  stright  line  trajectory). 

Let  a  plane  sounding  wave  (p0=0,  zQ-*  oo)  strike  the  object. 
Then,  by  calculating  the  integral  (11)  at  the  of  high  k  (kr^/2C»1 ) 
by  the  saddlepoint  method  (assuming  that  V(s,w)  =  -4tuC®1  in 
accordance  with  smooth  perturbation  method  approximation),  we 
obtain 

qz(p,k)  =  Jq(r,k)dz  -  (k/27cC)2  V(pR)  (2£/k)  (-1%)=  2ikffi1(pR). 

(12) 

To  demonstrate,  the  advance  of  the  complex  phase  ®  at  the  object 
in  the  geometrical  optics  approximation  with  a  linear  ray  is  equal 
to  the  integral  of  complex  refractive  index  n(r,co)  =  v  1-q/k2 

$  =  $0  +  $  +  ...  =  ik  J  n(r,o))  dz  = 

=  ik  J  v  1-q/k2  dz  ^  ik  J  (l-q/2112)  dz 

From  this,  for  the  first  approximation  of  the  complex  phase  of  the 
field  we  have  ^  “  -iq^/2k,  which  agrees  with  (12),  i.e.,  at  the 
high-frequency  limit  linear  integrals  of  type  (11)  are  only 
dependent  on  the  field  on  the  ray.  This  means  that  the 
approximation  used  here  also  contains  the  geometrical  optics 
approximation  with  linear  trajections  at  the  high-frequency  limit. 
The  result  of  the  limitation  on  s  of  the  region  where  the  field  is 
recorded  is  that  the  function  qz(p,w)  reconstructed  from  limited 


data  on  ©1  (s)  will  be  smoothed.  The  null  space  of  the  equation  is 
determined  by  the  parameters  of  the  recording  system  and  defines 
the  Rayleigh  resolution  limit  m;  then  qz(p,o))  is  an  integral 
over  a  ray  of  finite  "thickness". 

To  reconstruct  the  2D  structure  of  ionospheric  irregularities 
according  to  the  transformation  (11)  we  need  the  information  about 
the  approximate  coordinate  of  the  scatterer  zB.  Let  us  discuss 
briefly  a  method  for  determination  of  a  distance  to  the  scatterer, 
which  is  essential  for  irregularities  reconstruction  by  integral 
transformations.  If  only  an  approximate  value  of  the  distance  is 
known  then  there  appear  errors  in  s  and  C,  and  hence,  the 
reconstructed  function  is  destorted.  Let  us  consider  an  example 
where  a  source  Is  placed  on  board  of  a  satellite  and  the 
corresponding  distortions  arising  In  this  case.  The  distance 
between  the  satellite  and  the  receivers  (Zq  -  zR)  is  known  quite 
precisely,  but  in  determination  the  distances  "source  -  scatterer" 
(Zq  -  zg)  and  "scatterer  -  receivers"  (zg  -  zR)  a  systematic  error 
A  is  present:  zB  =  zg  —  A.  It  may  be  shown  that  reconstruction 
errors  are  small  for  such  values  of  A,  which  correspond  to  the 
longitudinal  resolution  of  the  measuring  system.  This  conclusion 
is  also  confirmed  by  numerical  simulations  [13. 

It  easy  to  show  that  when  the  error  A  is  less  than  the  longi¬ 
tudinal  resolution,  A  «X  2£2/(max  s)2,  the  reconstructed 
function  is 

k 

qz(p)  “  qz(x/v  y/V  exp(12 (y2~  x2))>  (13) 

where  ux  =  1  -  A/(zg-zR),  u^  =  1  +  A/(zQ-zs).  This  equation  makes 
it  possible  to  determine  the  distance  to  the  scatterer  by 
performing  reconstruction  with  different  values  of  A.  At  the  true 
location  of  the  scatterer  (A=Q)  the  "phase  front"  curvature  of  the 


function  under  reconstruction  changes  to  the  opposite  with 
respect  to  each  coordinate.  These  problems  are  considered  in 
detail  in  [13  and  partly  in  [63-  There  one  can  find  the 
generalisation  of  the  method  for  the  ionosphere  with  the 
stratified  background;  the  method  for  determination  of 
inhomogeneities  parameters  ("mass",  space  coordinates,  size  and 
other  moments  of  potential  q)  by  means  of  a  small  number  of 
receivers  and  solutions  of  IP  due  to  data  of  back-  scattering. 

It  is  possible  to  consider  the  more  general  case,  namely,  the 
strong  scattering  case.  For  this  case  in  accordance  with  (7,10)  we 
reconstruct  the  projection  of  the  product  of  the  potential  and 
the  normalized  field  Q=qU,  namely, 

J  qU  dz  =  Qa(p)  <14> 

It  should  be  stressed  that  in  general  the  normalized  field  U[q3 
depends  on  the  potential  and  on  the  direction  of  the  sounding  wave 
incident.  So  it  is  impossible  to  solve  the  problem  of  the  product 
qU  reconstruction  by  means  of  the  direct  transformation  of  a  set 
of  projections  (14).  It  is  possible  to  suggest  iteration 
procedures  of  the  ISP  solution  to  the  equation  (14)  directly.  But, 
if  the  scale  of  ionospheric  inhomogeneities  more  than  sounding 
wavelength  in  103-  104  times,  it  is  even  complicated  to  solve  the 
direct  scattering  problem.  So  for  such  cases  it  is  necessary  to 
use  asymptotic  methods.  The  special  asymptotic  methods  of  the 
solution  of  the  ISP  for  the  case  of  strong  and  large  scatterers 
were  developed  in  [1,183.  For  example  the  following  asymptotic 
representation  of  normalized  field  was  obtained 

i  s  1  r  d2q(p+z)  _3 

U(p,z)  =  exp[-  -  fq(p+z’)ds  +  — ~  (z-z*)— ^5-  dz'  +  0(k  J)1 

2k  -co  4k  l  dp1- 


(15) 


This  formula  is  correct  for  infinitly  smooth  potential.  The 
corresponding  additional  terms  in  the  sum  (15)  will  appear  in  the 
presence  of  derivatives  discontinuities.  The  formula  (15)  permits 
us  to  construct  the  simple  iteration  procedure  of  the  ISP 
solution,  when  it  is  not  necessary  to  solve  numerically  neither 
the  direct  scattering  problem,  or  the  ISP.  The  given  first  guess 
of  the  potential  (a  priori  information  or  transformation  under 
the  condition  of  weak-scattering)  defines  the  approximation  of  the 
normalized  field  functional  U[q1  ]  according  to  (15).  Using  this 
approach  it  is  possible  to  solve  the  usual  tomographic  problem  of 
the  reconstruction  q(r)  (14)  with  known  "weight"  U  due  to  Qa(p) 
under  the  different  directions  of  the  sounding  wave  incident. 
Later  on  the  iteration  procedure  with  the  potential  approximations 
obtained  is  repeated. 

The  formula  representation  (15)  permits  us  to  formulate  the 
conditions  of  applicatability  of  the  weak-scattering 
approximations.  According  to  this  condition  the  first  term  (^k-1 ) 
and  all  successive  terms  of  the  expansion  (15)  have  to  be  small. 
The  restriction  on  the  first  term  gives  us  the  well-known 
condition  of  the  BA  applicatability. 

3rm/2k  «  1  (16) 

Here  q  is  the  typical  (mean)  value  of  the  potential,  rm  is  the 
maximum  size  of  the  inhomogeneity.  In  the  case  when  the  wavelength 
of  the  sounding  radiation  \=2m  and  the  value  of  disturbances  of 
the  ionosphere  electron  density  N  ~  1 01 1  el/m3  the  condition  (16) 
leads  to  the  inequality  rm«  2  km.  The  disturbances  of  the 
electron  density,  scattering  radiowaves,  have  to  be  measured 
relative  to  the  background  of  the  regular  stratified  ionosphere. 
So,  if  the  background  value  N  «  1 Q1 1  el/m3  and  the  disturbances  are 


sufficiently  strong  (10%),  the  BA  is  valid  up  to  the  sizes  of 
dozens  tan.  But  for  the  case  of  the  main  ionosphere  maximum  (N  ™ 
10  m  )  the  scale  restriction  is  significant  and  it  is  necessary 
to  take  into  consideration  the  first  term  of  (15).  The  second  term 
of  (15)  for  simple  scatterers  without  internal  structure  is  less 
(kr^  »  1)  than  the  first  one.  In  the  case  when  there  are  some 
internal  details  inside  the  scatterer  with  the  scale  a»X  and  the 
typical  value  of  the  potential  q',  the  second  term  can  be  more 
than  the  first.  Hence  the  other  necessary  condition  of  the  BA  is 

the  condition  that  the  second  term  with  the  transverce  derivative 
is  small 

q’rm/(4k2a)  <c  1.  (17) 

In  many  practical  cases  it  is  sufficient  to  take  into  account 
the  first  term  of  (15)  only,  but  the  scatterer  is  not  the  Bom's 
scatterer  according  to  (16)  and  the  inequality  (17)  is  valid.  Such 
conditions  allows  us  to  get  the  simple  analytic  formula  [1,14], 
which  connects  the  function  Qs  reconstructed  due  to  the 
experimental  data  and  the  projection  of  the  scattering  potential 

qs  (is). 

Qa(p)  =  2ik  [exp(-iqs(p)/2k)  -  1]  (18) 

This  formula  is  the  basis  of  the  tomography  reconstruction  for  the 
most  of  ionospheric  applications.  Really,  the  inhomogeneities  with 
the  scale  about  ten  or  more  kilometers,  the  deviations  of  which 
are  about  10%  with  respect  to  the  background,  are  the  strong 
inhomogeneities  (the  left  hand  side  of  (16)  is  more  or 
approximately  equal  to  6,  when  X=  2  m),  but  the  inequality  (17)  is 
valid  for  the  internal  details  with  the  scale  more  than  hundreds 
of  meters  and  with  the  electron  density  disturbances  about  100% 
compared  with  the  background. 


The  uniqueness  of  solutions  of  the  similar  ISP  was  considered 
many  times  earlier  Cl -31.  Some  exact  results  are  known.  For  the 
problems  in  question  the  uniqueness  of  the  transformation  qz  (11) 
and  the  following  reconstruction  (14,18)  of  projections  can  he 
easily  proved  for  the  finite  volume. 

The  results  described  above  allow  us  to  get  2D  projections  of 
localised  ionospheric  inhomogeneities  with  the  help  of  the  moving 
satellite  transmitter  and  the  transverse  arrow  of  Earth  receivers. 
If  there  are  some  arrows  of  receivers  separated  by  distances  about 
hundreds  kilometers  along  the  line  of  the  satellite  flight,  it  is 
possible  to  obtain  some  2D  projections.  Having  a  set  of  2D 
projections  one  can  reconstruct  the  3D  structure  also.  Such 
experiment  (with  some  arrows  of  receivers)  is  rather  complicated 
and  it  was  not  realized  yet.  Nowadays  experiments  connected  with 
the  reconstruction  of  the  2D  localized  inhomogeneities  were 
carried  out  [1,5].  This  is  the  tomography  problem  also  because  it 
is  possible  to  reconstruct  the  2D  inhomogeneity  projection  or  the 
2D  cross-section  of  its  spectrum  after  the  corresponding 
transformations  by  means  of  a  set  of  ID  field  records  in  each 
receiver. 

We  will  consider  the  examples  of  numerical  calculations  by 
means  of  eleborated  programs. 

Using  the  first  numerical  example,  we  will  illustrate  the 
necessity  for  considering  diffraction  effects  in  RT  reconstruction 
of  objects  with  dimensions  rm  comparable  to  the  dimensions  of  the 
Fresnel  zone  v  A-C*  As  was  demonstrated  above  (12),  diffraction 
effects  are  insignificant  where  krJ^C  *  1  80)4  can  significantly 
distort  the  results  where  krjj/2C  <  1.  Fig.  2a  shows  qz(p), 
consisting  of  two  irregularities  -  two  Gaussians  (~  exp(-p2/r|) 
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with  parameters  (r^/v  AC  =  5,  r2/v  AC  -  1/3,  for  the  A  2m, 
C=200km,  2 i\,  =  6.3km,  2r2  =  420m).  If  we  reconstruct  the  complex 
phase  after  the  emission  passes  through  one  large  irregularity, 
then  according  to  (12)  Im  ft,  is  proportional  with  a  high  accuracy 
to  the  function  qz(p)  from  a  large  irregularity.  However,  the 
imaginary  part  of  the  complex  phase  ft,  (Pig.  2h)  from  a  pair  of 
irregularities  hears  little  resemblance  to  Pig.  2a  due  to 
diffraction  effects  in  small  irregularity.  Pig.  3  shows  the 
similar  picture  for  the  two  irregularities  with  parameters 
(r1  /v~JZ  =  5,  r2/vTC  =  1/8).  The  pair  of  irregularities  is 
reconstructed  very  accurately  in  the  case  where  the  Fresnel 
Inversion  (11)  is  used. 

For  reconstruction  of  the  structure  of  ionospheric 
irregularities  the  questions  about  sampling  data  and 
reconstruction  procedure  are  significant.  The  sounding  radio 
signals  are  always  recorded  by  point  receivers.  The  signal  from 
each  receiver  recording  the  emission  of  the  moving  transmitter  on 
the  satellite  is  also  sampled  by  an  electronic  circuit  with  a 
fixed  interrogation  frequency.  Moreover,  the  integral 
transformations  as  a  result  of  which  the  structure  of  the 
scattering  irregularities  is  reconstructed  can  be  realized 
numerically  on  a  computer  only  after  the  appropriate  digitization. 
Therefore,  both  the  form  of  the  data  recorded  and  the  nature  of 
the  numerical  processing  and  the  reconstruction  of  the  objects 
lead  to  the  necessity  for  performing  purely  discrete  operations. 
However,  in  theoretical  examination  of  the  solution  of  the  inverse 
scattering  problem  the  immediate  switch  to  discrete  foimulas  is 
hardly  advisable  in  the  majority  of  cases.  Usually,  the  discrete 
analogues  to  the  reconstruction  formulas  are  significantly  more 


cumbersome,  and  besides  the  discrete  form  hinders  analysis  of  the 
results  obtained.  On  the  other  hand,  the  switch  from  continual 
relations  to  their  discrete  analogues  does  not  cause  difficultes 
in  the  majority  of  cases  and  is  performed  in  a  unique  manner.  So, 
continual  formulas  can  be  considered  a  compact  notation  of  their 
discrete  analogues  for  further  practical  reconstruction. 

let  us  consider  the  switch  to  discrete  analysis  of  the 
reconstruction  procedure  and  the  questions  associated  with  this 
switch  using  the  transformations  (8),  (11).  It  should  be  noted 
that,  the  methods  for  digital  processing  of  signals  and  fields 
have  been  developed  rather  well  and  covered  in  the  literature 
[19,203 . 

First,  we  switch  in  (8,  (11)  to  the  dimensionless  variables  F 
=  p/v  A,£  ;  S  =  s/vTC  *  normalized  to  the  radius  of  the  Fresnal 
zone.  Intriducing  the  quantities 
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from  (11)  we  obtain  the  pair  of  integral  Fourier  transforms  for 
them  [1  ] 

F  (P)  =  X  FV(S)  e'2i%SP  d2S 

_  -*->  (20) 
FV(S)  =  f  Pq(P)  e+2l7cSP  d2P 

By  virtue  of  the  definition  Fq(P)  is  a  finite  function  due  to 
the  finiteness  of  qz(p) .  Let  the  dimensions  of  the  object  not 
exceed  the  limits  of  the  segment  t-P0x,  P0xJ  in  axis  x  and  [-P0y, 
P0y]  in  axis  y,  i.e.,  the  carrier  of  the  function  Fq(P)  is 
contained  in  a  fixed  rectangle.  Then,  Fy(S)  will  be  an  infinite 


and  analytical  function  with  a  finite  spectrum  F  (P).  According  to 
the  Kotel’nikov  theorem  FV(S)  can  be  represented  in  the  form  of  an 
infinite  series  in  terms  of  sample  functions  with  the  sampling 
intervals  4SX  =  P0x/2,  ASy  - 

PV(S  ,S  )  =  F  ( - ,  — 

y  WT  2F0x  2P0y 


)  sine  (2^^-^)  sinc(2icP0ySy-7any) 


Substituting  this  representation  in  (20)  and  integrating  with 
respect  to  Sx,  Sy,  we  obtain  the  relation 


(21) 


yp>  - 


“x  "V  -2ix[PASm  +  P  AS  ml 

—  >  e  1  (2Poz  2Poy i;1 


mx  “y  2P0x  2POy 


Using  the  concept  of  signals  which  are  really 

indistinguishable  at  the  level  e  [21],  it  is  advisable  to 

introduce  a  finite  function  which  is  really  indistinguishable  at  a 

fixed  level  from  FV(S).  let  the  carrier  of  this  function  be 

contained  in  the  rectangle  [-S0x,S0x]  x  C-S0y,S0y].  Then, 

according  to  the  Landau-Pollack  theorem  [21,22]  the  approximate 

■*> 

dimension  of  the  set  of  all  functions  F  (S)  finite  at  a  fixed 
level  with  a  carrier  in  the  rectangle  [-SQx,S0x]  x  [-SQy,S0y]  and 
with  a  finite  spectrum  supp  Fq(P)  €  [-PQx,  PQxl  *  t-P0y,  PQy]  is 
close  to  =  (2S0x*2P0x)  x  (2S0y*2P0y) .  In  other  words,  due 
to  the  finite  accuracy  in  measurement  of  the  signals  it  can 
assumed  that  both  the  signal  and  its  spectrum  are  finite.  The 
limitation  on  the  carrier  of  the  spectrum  leads  to  a  finite 
sampling  frequency  and 

Pq(P),  APX  =  1/2Ste,  APy  =  1/2Soy, 

here,  the  number  of  sampling  intervals  on  the  Cartesian  axes  is 

equal  to  =  Vp  Ny  ^o/^y  =  2W^r 

From  this,  assuming  Px  =  nxAPz,  Py  =  nyAPy,  it  Is  simple  to  switch 


from  (19),  (20)  to  a  pair  of  discrete  transformations  for  Fy  = 
V(W  H  K  “  pq 


^q^X*1^ 


1  Nx/2-1  Ny/2-1 

^x^y  mx=_Nx/2  my=_Ny/2 


Fv(lW 


-2iic( 
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TT  Ny  ’ 
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(22) 


VW  = 


Nx/2-1  Ny/2-1 

nx=_Nx/2  ny=-Ny/2 


V^’V 


mxnx  +  IW 
-2iiu(  N  N  ) 
a  x  y 


Later  on,  the  fast  Fourier  transformation  algorithms  are  used  in 
the  numerical  modeling  and  the  numbering  of  the  sums  over  m^,  n^ 
is  shifted  from  0  to  Nx  -  1 ,  and  likewise  for  ny  y. 

The  dimensions  of  the  reception  region  SQx,  SQy  and, 
consequently,  the  sampling  intervals  APX,  APy  also  determine  the 
resolutions.  To  demonstrate,  the  minimum  resolvable  interval  6X  In 
x  is  equal  to  6X=  rXC,  APX=  7TC/2S0x  =  =  A-C/2Sm.  A  limitation  of 
the  reception  region  always  limits  the  resolution  or  the  sampling 
interval  and  the  function  with  a  finite  spectrum.  Therefore, 
generally  speaking,  with  a  limited  reception  region  it  is  possible 
to  switch  immediately  from  (21)  to  the  finite,  discrete  Fourier 
transform  (22). 

The  switch  from  other  continual  reconstruction  formulas 
containing  the  Fourier  and  Fresnel  transforms  to  their  discrete 
forms  is  also  made  according  to  this  scheme. 

The  practical  realization  of  reconstruction  of  the  structure 
is  governed  to  a  great  extent  by  the  resistance  of  the 
reconstruction  procedure  to  the  various  types  of  noises  and 
distortions  which  unavoidably  appear  in  measurements  and 
processing  of  experimental  data.  Let  us  illustrate  the  influence 


of  noises  and  distortions  on  reconstruction  of  the  structure  on 
the  examples  of  numerical  simulation.  Briefly  discribe  a  modeling 
scheme  using  the  example  of  two-dimensional  Fresnel 
reconstruction. 

The  modeling  scheme 

^  ^  <v  ^ 

qz(p)  =>  P  (P)  -  PT(S)  *  V(S)  =■  .....  -  qz(p) 

^  I 

I 

I 

I 

noises,  distortions 

The  function  q  (p)  characterizing  the  two-dimensional  structure  of 

*2i  r 

the  electron  concentration  irregularities  and  the  effective 
collision  frequency  is  defined.  Then,  the  quantity  Fq(P)  is 
calculated  using  formula  (19).  After  performing  the  discrete 
Fourier  transform  (22)  corresponding  to  (20),  Fy(S)  is  obtained 
and,  thus,  the  field  V(S)  in  a  discrete  grid.  This  stage  of 
modeling  is  related  to  the  direct  scattering  problem,  where  it  is 
assumed  that  the  scattering  is  calculated  well  within  the 
framework  of  the  corresponding  approximations.  Further,  the  noise 
can  be  added  to  the  data  on  the  field  obtained,  after  which  the 
inverse  discrete  Fourier  transform  (22)  corresponding  (20)  is 
performed,  and  the  result  of  reconstruction  qz(p)  is  obtained  as  a 
result.  The  influence  of  distortions  is  modeled  in  asimilar 
manner:  in  the  inverse  transformation  stage,  the  parameters  of  the 
transformation  are  changed  or  ’’distorted". 

The  results  of  numerical  modeling  of  the  influence  of  noises 
have  demonstrated  the  stability  of  the  reconstruction  procedure. 
Fig. 4  shows  two  actual  irregularities  of  qz(p)  -  two  Gaussians 
with  various  dimensions  on  a  64x64  grid.  The  total  size  of  the 


image  frame  in  units  of  the  scale  of  the  Fresnel  zone  V  A,C  is 
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equal  5x5.  After  calculating  the  complex  phase  of  the  field 
scattered  by  these  irregularities,  the  data  on  V  were  disturbed  by 
noise.  Then,  reconstruction  was  performed;  the  quality  of  it  was 
sufficient  even  with  comparatively  large  errors  in  the  measured 
field  [233.  For  example,  Fig. 5  shows  the  result  of  reconstruction 
of  these  irregularities  according  to  data  on  with  additive 
complex  Gaussian  noise  tiaving  a  variance  of  0.05  of  maximum 
amplitude  of  the  change  in  the  real  and  imaginary  parts  of  complex 
phase  of  the  field.  Fig. 6.  illustrates  reconstruction  with  double 
the  noise  (a  variance  level  of  0.1).  Another  numerical  example: 
actual  irregularities  -  "cylinder”  and  "parabola”  are  represented 
in  Fig.T.  The  result  of  reconstruction  of  this  model  structure 
with  additive  noise  (a  variance  level  of  0.05)  is  shown  in  Fig.8. 
Fig. 9  illustrates  the  reconstruction  with  double  noise  (0.1).  If 
we  perform  preliminary  processing  of  the  data,  then  even  higher 
noise  levels  will  not  noticeably  influence  the  results  of 
reconstructions . 

Numerical  estimation  of  the  level  of  influence  of  the  noises 
is  of  interest  for  practical  applications.  It  is  advisable  to 
characterize  the  level  of  Influence  of  the  noises  in  the  metric  C 
and  L2.  In  a  discrete  reconstruction  procedure,  it  is  better  to 
use  their  discrete,  normalized  analogues,  i.e.  we  will  estimate 

rj 

the  difference  between  the  function  f  reconstructed  from  noisy 
data  and  real  function  f  by  the  numbers 


*  maxlf.-f^ 

p.(f,f)  =  -1 -  ;  p  2(f,f) 

0  max  |£,  |  l2 

i 


(fi-f± r 

i  1  1  1/2 

[  -  3  ,  (23) 
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where  summation  and  selection  of  maximum  are  performed  for  all 
increments  of  recorded  function. 


Zi 


The  dependence  of  the  normalized  disturbance  of  the 
reconstructed  two-dimensional  structures  (i.e.,  the  difference 

rv 

between  the  true  real  qz(p)  and  reconstructed  qz  in  the  metrics 

p  and  p  2)  on  the  amplitude  A  of  normalized  (to  the  maximum 
1 

amplitude  of  the  change  in  the  field)  complex  Gaussian  noise  is 
shown  in  Tabl.1  for  the  irregularities  in  Fig.7.  It  is  clear  that 
the  normalized  deviations  of  the  reconstructed  two-dimensional 
structures  are  comparable  to  the  noise  level.  The  absolute 
disturbances  of  the  real  part,  although  the  true  function  is  real. 

Distortions  of  irregularities  reconstructed  by  radio 
tomographic  methods  appear  due  to  imprecise  determination  of  the 
distance  to  the  irregularities  and,  as  a  result,  incorrect 
knowledge  of  the  given  parameters  of  the  integral  transformations. 
The  Influence  of  distortions  was  modeled  by  the  scheme  described 
above.  The  model  in  Fig.  10  (two  irregularities:  "cos”  and  "cosA2" 
ellipse)  is  the  object  of  reconstruction.  The  imprecisely 
determined  coordinate  of  the  scatterer  zg  =  zg  -  A  with  an  error  A 
was  substituted  in  the  parameters  of  the  Fresnel  transformation. 
Numerical  values  zQ  -  zg  =  700  km,  zg  -  zR  =  300  km,  X  =  2  m 
typical  for  satellite  sounding  were  used  in  the  calculations.  The 
dimension  of  the  Fresnel  zone  here  7~XC  *  0*65  km,  the  frame 
dimension  was  6, 4x6, 4  km  and  the  transverse  resolution  Q  was  of 
the  order  of  a  hundred  meters.  Fig.  11 -15  provide  reconstructed 
two-dimensional  structures  with  various  values  of  the  error  in 
determination  of  the  coordinates  A  measured  in  kilometers;  A  is 
successively  equal  to  -10,  +20,  -20,  +30,  +40  km.  It  is  clear  that 
an  error  of  10  km  has  little  influence  on  the  reconstruction. 
Errors  A  of  the  order  of  20  km  distort  the  object,  but  the 
characteristic  details  of  the  object  are  reconstructed  fairly. 
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Only  errors  exceeding  30-40  Km  noticeably  distort  the 

reconstructed  irregularity  and  may  lead  to  incorrect  conclusions  on 
the  structure  of  the  object.  The  values  of  the  errors  in 
reconstruction  in  normalized  metrics  are  shown  in  Tabl.2.  The 
modeling  was  performed  for  two  Gaussians.  The  analogy 

dependence  of  error  reconstruction  on  the  A  is  shown  in  the  graphs 
in  [1].  The  reconstruction  errors  have  a  clearly  pronounced 
minimum  and  are  minor  («$  0,1  t  0,2)  within  ±  (5  -  10)km.  We  will 
note  that  dimension  of  this  minimum  agrees  in  order  of  magnitude 
with  the  longitudinal  resolution  8Z  of  the  system;  here,  the 
transverse  resolution  8  ~  A./6  ~  100  m,  the  aperture  angle  9  ~0,02 
h  8  ~  A./92  ™  5  Km.  Therefore,  determination  of  the  distance  to 

scatterer  with  accuracy  of  the  order  of  the  longitudinal 
resolution  of  the  recording  system  is  sufficient  for  high-quality 
reconstruction.  The  method  for  determining  the  distance  according 
to  the  change  in  the  curvature  of  the  "phase"  of  the  reconstructed 
function  maKes  it  possible  to  achieve  this  accuracy  in 
insignificant  noises  Cl  3 . 

On  reconstruction  of  qz  using  the  formula  (18)  the  main 
obstacle  is  that  Qqz  is  a  periodical  function  of  qz,  therefore,  a 
real  part  of  qz  is  not  unambigiously  defined  by  Qz,  which  we  get 
from  experimental  data.  Since  that,  making  use  the  formula  (18) 
one  cannot  immediately  reconstruct  scatterers  that  exceed  4iuk.  It 
should  be  noted  here  that  it’s  much  more  than  BA-based 
reconstruction  procedures  are  applicable  to.  Namely,  we  have  found 
using  numerical  simulations  that  BA  holds  for  scatterers 
satisfying  the  inequality 


qz/2k  <0.1. 


In  principle,  one  can  reconstruct  more  intensive  scatterers, 
assuming  that  Qz  is  zero  if  qz  crosses  iuk  level,  or  if  any 
additional  information  about  the  scaterer  structure  is  available. 

We  have  derived  an  algorithmic  procedure  for  finite 
scatterers,  which  allows  one  to  reconstruct  scatterers  exceeding 
4iuk.  We  have  seen  that  the  main  obstacle  which  restricts  this 
procedure  capabilities  is  the  finite  number  of  receivers  which 
register  the  scattered  field.  From  the  formal  mathematical  point 
of  view,  that  simply  means  we  get  only  finite  number  of  Fourier 
expansion  terms  ,  i.e.  we  deal  with  the  truncated  Fourier  series 
of  Qz. 

Thus  to  obtain  qz  we  have  to  derive  it  from  a  truncated 
Fourier  expansion  of  Q  ,  which  in  turn  is  a  nonlinear  transcendent 
function  of  q  .  This  makes  us  use  the  numerical  simulations  to 
investigate  the  domain  of  applicability  of  this  procedure,  because 
immediate  analysis  of  that  algorithm  is  extremely  difficult. 

We  have  performed  such  an  investigation  of  this  procedure 
applying  it  to  various  scatterers.  The  examples  of  such  numerical 
simulations  of  reconstruction  of  ID  projections  of  2D  scatterers 
are  presented  on  Fig. 16-1 8.  The  scatterer  -  receiver  distance  was 
assumed  to  be  300  km,  the  sounding  wavelength  was  2m.  Thick  curves 
correspond  to  initial  projections  qz  which  had  to  be 
reconstructed,  thin  curves  are  the  reconstructed  images  of  qz 
using  the  procedure  described  above.  Just  for  checking,  the  Qz 
patterns  are  also  shown  (small  amplitude  fast  oscillating  curves 
of  normal  thickness).  On  the  Figure  16,  the  scatterer  is  not  too 
strong,  so  two  curves  simply  coincide  (  the  projection  is  exactly 
reconstructed  ).  On  the  Fig.  17-1 8,  where  stronger  patterns  are 
given,  one  can  see  the  discrepancy  between  exact  image  and 


reconstruction.  One  can  see  very  easy  that  mistakes  in  the 
procedure  work  arise  mostly  on  those  fragments  of  scatterers  where 
the  transverse  derivative  idqz/dxi  =  !dQz/dxi  is  big  enough,  so 
the  oscillations  of  Qz  become  so  fast  that  Qz  is  no  longer 
approximated  by  the  truncated  Fourier  expansion  well  enough.  They 
may  cause  the  loss  of  zeros  of  Qz,  which  correspond  to  increments 
or  decrements  of  4itk.  The  numerical  simulations  using  wide  variety 
of  scatterers  have  confirmed  this  conclusion.  Moreover,  they  have 
shown  that  only  the  transverse  derivative  absolute  value  idqz/dx! 
is  critical  for  the  correct  reconstruction.  Of  course,  it  causes 
the  limitation  on  the  maximal  intensity  of  the  scatterer,  since 
for  continuous  function  qz(x)  !  qz  I  <  a  idqz/dxi.  For  this  reason, 
the  problem  to  calculate  this  boundary  value  for  idqz/dxi  had  also 
been  posed.  The  value  which  we  get  estimating  the  supremum  of  the 
truncated  Fourier  expansion  derivative  (it  is  bounded  as  a 
derivative  of  polynom)  is  much  more  than  actual  value.  Fig.  19 
shows  the  normalized  derivative  value  times  2icN  plotted  versus  N, 
where  N  is  the  number  of  terms  retained  in  the  Fourier  series.  We 
assumed  that  well-approximating  reconstruction  it  that  one  which 
has  the  maximal  normalized  discrepancy  between  exact  and 
reconstructed  images  less  than  0.1. 


The  description  of  Fortran  Programs 
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System  Requirements 


-  Computer: 

-  Operating  System: 

-  Memory: 


-  Hard  Disk  Space: 

-  Software: 


IBM  AT  or  compatible  (with  coprocessor) 
MS-DOS  or  PC-DOS  version  3.0  and  later 
at  least  Extended  memory  4  Mbytes 
(depends  on  geometry  and  type  of 
approximation  reconstructed  function) 

4  Mbytes 

compiler  1.4e  and  linker  2. 2d 
NDP-F0RTRAN-386(c)  MicroWay  or  later 


1.  Program  <DT_DIR.for> 

This  program  solves  direct  problem,  namely,  determines  complex 
field  V(S)  by  means  of  Past  Fourier  Transformation  (FFT)  based  on 
the  model  structure  of  irregularities. 

Input  parameters  and  files: 

(the  determination  of  the  irregularities  parameters) 
file  <mod_par.tsk>  contains  the  following  parameters: 

«  line  1  in  file 


'wave  length  (km):' 
line  2 

0.002 

(DIV  in  program) 

'distance  to  satellite  (km):’ 
line  3 

1000. 

(ZZO  in  program) 

'distance  to  irregularities  (km):' 

300. » 

(Z  in  program) 

These  parameters  you  can  change  in  the  file. 

NY  -  number  of  receivers 

MNY  -  the  degree  of  2  for  NY  (namely,  2AMNY=NY) 

NX  -  number  of  discrete  intervals  on  the  axis  X 
MNX  -  the  degree  of  2  for  NX  (namely,  2AMNX=NX) 

RMN  -  the  radius  of  the  Fresnel  zone. 

PX.PY  -  the  limits  of  segments  t-PX/2,  PX/2]  in  axis  X  and 
[-PY/2,  PY/2]  in  axis  Y;  object  (model  structure  not 
exceed  these  limits;  variable  P) 

SX.SY  -  the  limits  of  rectangle  [-SX/2,  SX/2]x[-SY/2,SY/2] 
(variable  3) 


Si 

YR  -  the  size  of  the  receiving  system  (km) 

X02  -  the  size  of  aperture  ’’synthesized.’*  by  satellte  (km) 

DXX,  DYY  -  the  size  of  frame  (km*km) 

NMOD  -  number  of  model  structure  (1-10): 

1  -  the  gaussians 

2  -  ’’homogeneous”  ellipse 

3  -  ’’parabolic”  ellipse 

4  -  ”  2  +  3  ” 

5  -  ”  2  ’’  with  absorption  10% 

6  -  ”  3  ”  with  absorption  20% 

7  -  ”cos”  ellipse 

8  -  "cosA2”  ellipse 

9  -  ”  7  +  8  ” 

10  -  arbitrary 

Nmod  is  introduced  from  the  screen 
Const5=0.1  (namely, 10%)  -  constant  for  NM0D=5 
Const6=0.2  (20%)  -  constant  for  NM0D=6 

Z1 ,Z2,Z3,UG  - 

A1 ,B1 ,A2,B2,X0,Y0,Y1 ,A3,FPC,FPS  parameters  for 

AEEP.BEEP 

APAR.BPAR  models 

AELP1 ,BELP1 ,APAR1 ,BPAR1  - 

RN  -  complex  field 

Subroutine  <SMAXMI>  calculates  min  and  max  of  array 
Subroutine  <FFT>  calculates  1-D  Fast  Fourier  Transformation 

Output  files: 

<M0DEL.GRD>  -  contains  the  model  structure  of  irregularities  and 
the  number  of  the  model  (NMOD)  and  constans  VMD  and 
ZRD  (for  determination  of  reconstruction’s  errors) 
<V_real.grd>  -  the  array  of  real  parts  of  the  field  RN 
<V_imag.grd>  -  the  array  of  imaginary  parts  of  the  field  RN 
Execution 
f77  dt_dir.for 
RUN 

dt_dir.exp 


2.  Program  <DT_INV.for> 

This  program  solves  inverse  scattering  problem  for  weak- 
scattering  irregularities  in  the  case  when  noise  is  present  in  the 
measurements  of  the  field  and  calculates  the  errors  of  the 
reconstruction  in  metric  C  and  12  in  dependence  on  level  of  noise. 

Input  parameters  and  files: 

Parameters  <  NX,  MNX,  NY,  MNY,  PX,  PY,  SX,  SY  >are  similar  to  same 
parameters  of  program  <DT_DIR.P0R>. 

FILES  (from  program  <DT_DIR.for>): 

<M0DEL.GRD>  -  input  file  of  model  structure. 

<V_real.grd>  -  the  array  of  real  parts  of  the  field  RN. 
<V_imag.grd>  -  the  array  of  imaginary  parts  of  the  field  RN. 
Parameter  <Error>  -  the  level  of  noise  (0  <  error  <  1 ) 

Error  is  introduced  from  the  screen. 

ConstS,  Const6  -  are  similar  to  same  parameters  of  program 
<DTJDIR.for> 

Subroutine  <SMAXMI>  calculates  min  and  max  of  array 
Subroutine  <FFT>  calculates  1-D  Fast  Fourier  Transformation 
Function  <RAN>  -  determinates  random  values  in  [0,1]. 

K33  -  the  constant  for  function  <RAN> 

Output  files: 

REC_er.GRD  -  the  reconstruction 

er_rler.dat  -  contains  the  errors  of  reconstruction  in  metric  C 
and  12 

Same  errors  are  shown  on  the  screen. 

Execution 
f77  dt_inv.for 

RUN 

dt_inv.exp 

3.  Program  <DT_ERR.for> 

This  program  solves  inverse  scattering  problem  for  weak- 
scattering  irregularities  in  the  case  when  distortions  are  present 
in  the  measurements  of  the  field  and  calculates  the  errors  of  the 
reconstruction  in  metric  C  and  12  in  dependence  on  distortions. 
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Input  parameters  and  files: 

Parameters  <  NX,  MNX,  NY,  MNY,  PX,  PY,  SX,  SY  >  are  similar  to  same 
parameters  of  program  <DT_DIR.F0R>. 

FILES  (from  program  <DT_DIR.for>): 

<M0DEL.GRD>  -  input  file  of  model  structure. 

<V_real.grd>  -  the  array  of  real  parts  of  the  field  RN. 
<V_imag.grd>  -  the  array  of  imaginary  parts  of  the  field  RN. 
<mod_par.tsk>  is  similar  the  same  file  in  program  <DT_DIR.for>. 
Parameter  <DEL>  -  the  error  of  Irregularities  coordinate  in  km 
DEL  is  introduced  from  the  screen. 

Subroutine  <FFT>  calculates  1-D  Fast  Fourier  Transformation 

Output  files: 

REC_cr.GRD  -  the  reconstruction 

er_rlcr.dat  -  contains  the  errors  of  reconstruction  in  metric  C 
and  12 

Same  errors  are  shown  on  the  screen. 

Execution 

f77  dt _ err. for 

RUN 

dt_err.exp 


4.  The  solution  of  the  statistical  radiotomography  problem 
(statistical  inverse  scattering  problem) 

Rather  frequently,  the  ionosphere  contains  entire  regions 
filled  with  large  numbers  of  irregularities  of  various  dimensions. 
Such  a  state  is  typical  for  the  equatorial  and  polar  regions, 
especially  at  night.  In  this  case  it  is  advisable  not  to 
reconstruct  individual  realization  of  such  disturbed  ionosphere, 
but  rather  to  pose  the  problem  of  reconstruction  of  the 
statistical  characteristics  of  the  turbulent  random  ionosphere 
such  as  correlation  function  or  spectrum  of  fluctuations.  So,  such 
reconstruction  of  random  ionosphere  statistical  properties  by  the 
measured  field  statistics  it  is  reasonable  to  name  the  statistical 
ISP  or  the  statistical  RT  problem  [24,25]. 

We  consider  the  reconstruction  of  ionosphere  electron  density 
fluctuations  spectrum  by  means  of  satellite  radioprobing  data. 
Integral  equations,  which  relate  the  measured  field  to  the  medium 
structure  are  probably  the  most  convenient  and  adequate 
mathematical  technique  for  the  tomography.  In  this  paper  we  shall 
show  how  to  derive  such  equations.  For  simplicity  it  is  reasonable 
to  use  two  frames  of  reference.  The  first  "global”  frame  of 
reference  r=(x,y,z)=(p,z),  its  origin  is  related  to  the  receiving 
system  as  in  the  previous  section  (Fig.1).  The  origin  of  the 
second  frame  of  reference  'R=(X,y,Z)=(P,Z)  is  reasonable  to  locate 
at  one  of  the  transmitter  positions  and  the  Z  axis  to  direct  to 
the  center  of  the  receiving  system.  As  the  transmitter  moves  it  Is 
convenient  to  use  several  such  "local"  frames  with  different 
orientations.  Such  frames  are  introduced  to  make  derivations 
shorter,  because  the  probing  wave  scattering  in  each  local  frame 
may  be  treated  as  "almost  forward"  small-angle  scattering.  The 
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lower  boundary  of  the  layer  of  irregularities  is  z^,  the  upper 
zu  in  the  global  frame  r.  After  performing  the  Fresnel  expansions 
of  exponents  in  the  (5)  for  small-angle  scattering,  we  obtain  the 
integral  equation  for  U  (7).  This  equation  is  written  in  the  local 
reference  frame,  S(Z-ZQ)  =  (Z'-Zq)H  +  (Z-Z')ft  ;  Z>(Z-Zq)  = 
(Z’-ZQ)  (Z-Z' ).  For  making  formulas  shorter  we  assume,  that  the 
transmitter  is  located  exactly  at  the  local  frame  origin,  i.e. 

=  (P0,Z0)=(  0,0,0). 

To  reduce  calculations  and  to  use  known  results  we  transform 

-► 

(7)  to  the  parabolic  equation  in  new  variables  £  =  1/2;,  o  =  P/Z, 
!'=  i/Z’.  Then  (S-P' )2/ V  =  (o-l’P'  )2/ and  U  satisfies  the 
differential  equation 

(~2ik  —  +  A  -  r2q)  U(£,o)  =  0  (24) 

dt,  ° 

Equation  (24)  is  derived  from  (7)  by  differentiating  and  making 
use  of  the  following  relation  for  the  fundamental  solution  of  the 
Schrodinger  operator  [24].  Having  derived  the  parabolic  equation 
we  use  the  methods  developed  to  obtain  equations  for  the  first  and 
the  second  moments  of  U  [26,27].  Let  us  assume  also  that  the 
random  field  q  to  be  Gaussian  and  5-correlated  along  Z, 

B (7?1  ,P2)=  <  q (P1  ,Z1  )  q(P2,Z2)>  =  6(Z2~Z1  )  kQ{V2-V1  )  (25) 


In  the  similar  way  in  the  variables  (o,£)  it  is  possible  to  derive 
from  (24)  equations  for  the  second  moments  of  the  normalized  field 
7(7?)  =  E(7?)/<E(7?)>, 

■4  -4 

r1  0(P,Z)  =  <7(P2,Z)  7*P1,Z)>  -  the  first  coherence  function, 

-4  *4  4 

T0  Q(P,Z)  =  <7(P2,Z)  7(P1,Z)>  -  the  second  coherence  function  of 
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the  second  order.  As  the  calculations  are  described  in  the  papers 

mentioned  above,  we  present  only  the  final  result  in  the  variables 

P+  =  (P^+P^)/2t  P  =  P2.-P1  and  for  statistically  homogeneous  layer, 

when  the  dependence  on  P+  is  absent. 

-  (26) 

Q  V  d  i 

(tk( —  +  -  — )  +  +  —  a  (P)i  r2>0  =  o,  r2>0(£,zu)  =  i 

dZ  Z  dV  4k 


The  integral  equation,  corresponding  to  (21 ) ,  has  the  form 


An(^) 


-  i  r  o  V  >  ^  ik  .  .  2 

r?  n(P,Z)  =  1  + -  d3ft’ - r2  Q {P',Z')  exp(7 > 

16iuk  J  V  ^’u  4Z/ 

The  following  formula  for  r1  -j  is  well  known  [26,271. 


(27) 


TU^P,Z)  =  exp(^  J 

Eguation  (27)  follows  from  (26).  The  resulting  integral  equation 
is  valid  also  at  PQ  ^  0,  ZQ  /  0. 

The  integral  equation  (27)  is  the  basis  of  the  statistical 
RT.  According  to  (20)  Aq(7*)  is  the  correlation  function 

projection,  /  B (P,Z)  d Z  =  A^(P).  Measured  wave  coherence  functions 
allow  us  to  determine  projections  of  the  complex  potential  q 
correlation  function  by  the  (27,28).  To  begin  with,  let’s  assume, 
that  irregularities  occupy  a  sufficiently  extensive  layer,  oblique 
to  the  probing  wave  propagation  direction.  To  reconstruct 
projection  A ^(P)  it’s  necessary  to  determine  the  layer  coordinate 
Za.  It  may  be  done  by  a  single  receiver.  The  special  procedure  for 

o 

determination  of  the  scattering  layer  ^-coordinate  by  equation 
(22)  was  developed  [1,241. 

In  many  situations  randomly  inhomogeneous  media  represent 


Z ’ 

Aq(p-2_) 


dZ’) 


(28) 


statistically  quasi-homogeneous  fields  with  slowly  varying 
statistical  properties,  i.e.  such  fields,  for  which  the 

correlation  radius  of  the  differential  argument  is  essentially 
less  than  the  scale  of  the  variance  a.  In  other  words,  the 
correlation  function  is  represented  as  a  product  B  {ft^  ,ft^ )  = 
a2  (ft)1C{hft,ft) ,  where  ft  -  (ft^  +  ft^  )/2,  Aft  =  ft^ -ftg '  ^ 
inhomogeneous  field  with  the  costant  correlation  coefficient  £(Ar) 
hut  with  varying  fluctuations  variance  a2  (ft)  may  he  called  an 
"additive"  field,  because  the  electron  density  changings  influence 
only  the  fluctuations  intensity,  but  the  correlation  coefficient 
dependence  on  the  summarize  argument  ft  is  absent.  Assuming  the 
field  q  to  be  6-correlated,  as  was  done  above,  it  is  possible  to 
obtain  similar  equations  by  the  same  technique.  The  correlation 
function  projection  Aq(P)  in  the  integrands  should  be  replaced  by 
o2(Z)  K(P).  Here  K(P)  =  X  JC(Aft)  d(AZ’)  is  the  correlation  function 
projection.  Then  formula  for  I\,  ^  and  the  equation  for  r2>0  have 
the  form 


r2,0(^Z)=  1  + 


—  J 

167ck  J 


d3ft’ 


o2(^,)K(^') 


k 


V 


r2j0(7?’)  exp(i W(S-P')Z) 


(29) 


•j  Z  2' 

r1  1  {P,Z)  =  exp(  -^2  J  o2(2’’)  K(p  2~)  )  <30) 

Zu 

The  equations  (29,30)  contain  the  unknown  function  a2  (ft).  So,  the 
problem  of  the  correlation  function  reconstruction  is  divided  into 
two  ones:  the  problem  of  reconstruction  of  the  fluctuations 
variance  o2 (ft) ,  and  the  problem  of  the  reconstruction  given 

the  measured  projections  K(£). 

The  formula  (30)  allows  us  to  determine  integrals  of  o2(2) 
along  the  direction  (Z  axis)  of  the  probing  wave  propagation.  If  a 
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set  of  different  integrals  is  given  for  a  region  filled  with 
irregularities  we  arrive  at  a  tomography  problem.  Usually  in  the 
experiments  there  exist  a  natural  limit  to  the  number  of 
projections  and  angular  range  of  probing  waves  propagation 
directions.  Hence,  the  problem  when  data  are  Known  in  the  small 
angle  range  should  be  solved.  A  similar  problem  was  solved  in  the 
ray  RT.  Information  about  the  fluctuations  intensity  a2  (ft) 
distribution  is  sufficient  to  reconstruct  a  set  of  correlation 
function  projections  given  by  measured  coherence  functions. 
Methods  of  the  statistical  RT  problem  solution  are  similar  to  the 
discussed  above,  with  the  substitution  of  the  product  c^WKtft) 
of  the  known  function  cr( li)  and  the  correlation  function 
projection  K(ft)  instead  of  A^.  A  set  of  projections  K(ft)  allows  us 
to  perform  tomography  reconstruction  of  the  correlation 
coefficient  £(A ft)  or  its  spectrum. 

With  the  help  of  examples  we  want  to  illustrate  the  results 
of  programs  and  numerical  pattering.  Reconstruction  of  electron 
concentration  fluctuations  spectra  was  made  in  two  steps: 
pattering  of  environment  fluctuation  spectra  and  recovery  of  the 
power  angle  spectra  of  the  dispersive  environment. 

I.  Pattering  of  environment  fluctuation  spectra.  In  the  given 
program  <ST-DIR.F0R>  the  account  of  accidental  realisations  of 
electron  concentration  power  fluctuation  spectra  is  made.  The 
pattering  scheme  is  given  of  tabi.3. 

1 .  The  sum  of  three  (of  five )  Gauss  lines  with  def f erent  width 
and  intensivity  was  taken  as  a  pattern  of  spectra  with  a  few 
maximums.  The  location  of  local  maxlmums  and  their  width  depend  on 
the  dimentions  of  the  area  (32*32, 64*64, et c. )  and  were  given  by 
the  quantity  of  accounts  from  the  beginning  of  the  area. The 
intensivity  of  local  maximums  was  given  in  the  percent  relativity 
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from  the  intensivity  of  the  center  maximum. 

2.  The  computation  of  correlation  function  projections  of 
dispersed  field  is  made  In  polar  coordinates  for  the  given 
reception  angles  In  condition  with  equation  (28).  The  projections 
of  environment  correlation  functions  RCp.cJ^)  calculated  with  the 
help  of  fast  Fourier  transformation  from  the  spectra  ^tp.c^)  for 
every  given  reception  angle  in  area  [0,i u/2]. 

3.  The  form  of  casual  realisations  of  the  dispersed  field  is 
made  in  two  steps.  At  first,  with  the  help  of  linear  congruent 
method  [28,29]  statlsticaly  independent  one-dimentional 
realisations  of  casual  numbers  with  proportional  spectra  were 
obtained.  Then,  with  the  help  of  spectra  method  of  reorganization 
[30]  the  prosess  with  given  correlation  function  r^Cp,^)  was 
formed. 

4.  The  selected  (casual)  correlation  functions  of  the  field 

rv 

l\,i  (p, 4^}  were  calculated  from  alone  realisations  of  the  field. 

5.  The  selected  projections  of  correlation  functions  of  the 

rv 

environment  R(p,<J>n)  were  calculated  by  logariphmation  of  the 

rv  rv 

function  l\j  .j  (p,<f>n).  The  appropriate  R(p,4>n)  casual  power  spectra 

rv 

of  environment  fluctuations  ^(p,^)  are  defined  with  the  help  of 
standard  procedure  fast  Fourier  transformation.  The  given  model  of 
environment  spectra  (in  Decart  coordinates)  is  given  as  the 
illustration  on  fig. 20, a.  Area  dimentions  are  64*64,  width  of 
central  and  local  maximums  on  the  level  of  0.5  their  intensivity 
is  1/8  from  linear  dimentions  of  the  area.  The  local  maximum 
lntensivety  is  1-^=0. 3IC  (Ic  means  intensivity  of  central  maximum). 
From  the  fig  20, b  one  can  see  casual  spectra  of  environment 
fluctuations,  formed  from  one  field  realisation  for  every  given 
reception  angle.  It  is  easy  to  see,  that  for  more  calculations 
reconstraction  of  angle  spectra  )  the  same  casual  function  need 


smoothing,  i.e.  period  of  Its  oscilatlons  Is  equal  to  the  account 
step  [291.  On  flg.21  statistically  middled  (by  50(a)  ond  100(b) 
realisations)  angle  spectra  of  the  intensivity  of  environment 
fluctuations  ,  took  by  the  given  below  scheme,  are  given. 

II.  Recovery  of  the  power  angle  spectra  of  the  dispersive 
environment . 

In  program  <ST~INV.F0R>  the  account  of  reconstructions  of 
environment  power  angle  spectra  from  the  information  In  finit 
reception  angle  area.  For  solving  the  problem  of  recovery  the 
method  of  function  decomposition  into  Lejandr  polynoms  is  used. 
The  system  of  linear  equations,  made  by  this  method,  is  solving  by 
Gordon  method  [301.  On  the  last  step  of  the  account  the  linear 
interpolation  of  the  function  $(p,<)>n)  in  Decart  coordinates  is 
made. 

On  fig. 23  one  can  see  reconstructions  of  power  spectra  of 
plasma  fluctuations  calculated  from  one  realisation  of  dispersed 
field  in  each  direction  of  the  reciption  of  the  signal.  The 
account  was  made  for  three-racurs  reception  system  (fig.22,b).  It 
was  supposed  that  angle  areas  makes  0-20°, 35-55°, 65-85°  from 
vertical.  On  fig. 23  the  model  and  reconstruction  of  its  casual 
realisation  (the  width  of  central  and  local  maximums  is  1/8  and 
1/10  from  the  linear  dimentions  of  area  (64*64), intensivity 
I-j=0.3Ic)  is  shown.  On  fig.24  the  results  of  the  same  account  for 
the  dimentions  of  local  and  central  maximums  of  1/16  and  1/8  is 
shown.  Analising  the  results,  we  can  see  that  the  best 
reconstructions  are  obtained  for  the  models  of  sufficiently  narrow 
spectra  (fig.25-28).  In  this  case  (fig.28)  the  recovery  of 
realisation  of  pattern  from  5  local  maximums  (width  is  1/20, 
relativity  of  intensivity  is  1:0. 3:0. 2)  is  sufficiently  good.  On 
fig. 29  the  reconstruction  of  realisation  of  pattern  from  5  local 


maximums  (width  1/16,  relativity  of  intensivity  is  1:0.3:0.2)  is 
given. The  results  of  the  accounts  also  show  that  weak  intensive 
local  maximums  (1:0.1)  recover  sufficiently  good  (fig.  30).  But 
the  essential  part  is  played  hy  the  width  of  the  spectra. On  fig. 31 
the  reconstruction  of  model,  which  appropriate  to  the  experimental 
observations [24] .  It  is  clear,  that  recovery  of  the  spectra  from 
one  realisation  allow  to  underline  the  particularities  of  its 
structure.  Comparison  analisys  shows  that  in  spectra 
reconstructions  the  locations  of  local  maximums  coincide  with  the 
original  but  the  intensivities  are  differ.  'It  future  thes 
difference  can  be  overcomed  by  using  of  correlative  filters. 

The  quality  of  reconstruction  very  much  depend  on  the 
quantity  and  width  of  reception  areas.  On  fig.32,a  reconstruction 
of  model  from  fig.24,  made  by  two-racurs  scheme (fig.22),  is  given. 
Reception  angles  area  -  0^30°  and  50^80°. *The  same  reconstructions 
are  made  for  the  experimental  model  (fig.30).  It  is  clear,  that 
3-racurs  sceme  is  prefered.  On  fig. 33  are  given  spectra 
reconstructions  for  three  reception  areas:  0-10°, 40-50°, 70-80°, 
which  are  in  connection  with  fig. 32.  It  is  clear  that  to  give  more 
quality  to  the  reconstruction  the  width  of  reception  area  should 

have  optimum  value. 


The  description  of  Fortran  Programs 


System  Requirements 


-  Computer: 

-  Operating  System: 

-  Memory: 


-  Hard  Disk  Space: 

-  Software: 


1.  Program  <ST-DIR.for> 

This  program  solves  direct  problem,  namely, 
determines  the  model  structure  and  calculates 
the  random  realizations  of  angle  power  spectrum. 

Input  parameters  and  files: 

M  -  size  of  frame 
NM  -  size  of  half  frame 
JMN  -  number  of  rays  in  CO;  tu/2  3  range 
KR1  -  local  maximum  size 
KR  -  central  maximum  size 
N00  -  spectrum  coordinate  of  local  maximum 
OTT  -  level  of  local  maximum  power 
array  HH  -  model  of  angle  power  spectrum  in  polar 
coordinates 

array  SPAD  -  correlation  function  model  of  scattering  field 
Subroutine  SMAXMI  calculates  maximum  and  minimum  of 

one-dimension  function 

Subroutine  SMAXNI  calculates  maximum  and  minimum  of 

two-dimension  function 

Subroutine  SPGEN  calculates  balance  multipliers 
Subroutine  RANFIL  calculates  random  fields (complex) 
Subroutine  CORPAT  calculates  random  correlation  function 
Subroutine  F0RT1  calculates  one-dimension  Fast  Fourier 

Transformation 

Function  RANI  calculates  uniform  random  number 


IBM  AT  or  compatible  (with  coprocessor) 
MS-DOS  or  PC-DOS  version  3.0  and  later 
at  least  Extended  memory  4  Mbytes 
(depends  on  size  of  frame  and 
namber  of  rays) 

4  Mbytes 

compiler  1 .4e  and  linker  2. 2d 
NDP-F0RTRAN-386(c)  MicroWay  or  later 
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Output  files: 

<RANDSP.GRD>  -  contains  random  realization  of  angle 
power  spectrum,  size  (JMN:M) 

Execution 

f77  st-dir.for 

RUN 

st-dir.exp 

3.  Program  <ST-INV.for> 

This  program  solves  inverse  problem  for  statistical  RT, 
namely,  calculates  reconstructions  of  power  spectrum  of  electron 
density  in  Decart  coordinates  for  varios  namber  of 
receivers. 

Input  parameters  and  files: 

Parameters 
M  -  size  of  frame 
NM  -  half  size  of  frame 
NN  -  half  size  of  frame  +1 
ITER  -  number  of  iterations 
JPI  -  number  of  reseivers 
JMN  -  numer  of  rays  in  [0,iu/2]  range 
JT1 ,  JT2,  JT3,  JT4, JT5, JT6  -  boundaries  of  angle  ranges 

of  reception  in  degreep 

NU  -  number  of  rays  in  angle  ranges  of  reception 
NXX  -number  of  Lejandr  cofficient 

PILE  <RANDSP.GRD>  contains  realization  of  angle  power  spectrum 

(OUTPUT  file  of  progpam  ST-DIR.for) 

OUTPUT  FILE: 

<  RECONS. GRD>  contains  the  reconstraction  of  realizations 
angle  power  spectrum  in  Decart 
coordinates , size (M : M ) . 


Execution 

f77  st-inv.for 

RUN 

st-inv.exp 


References. 


1 .  V.E. Kunitsyn,  E.D. Tereshchenko,  Tomography  of  the  Ionosphere , 

Moscow,  Nauka,  1991  (in  Russian). 

2.  V.E. Kunitsyn,  "Inverse  Scattering  Problem  in  the  Layered  Medium", 

Moscow  University  Physics  Bulletin,  Vol.26,  No  6,  pp. 33-40,  1985. 

3.  V.E. Kunitsyn,  "Determination  of  the  structure  of  ionospheric 

irregularities".  Geomagnetism  and  Aeronomy,  Vol.26,  pp. 75-81 ,  1986 

4.  V.E. Kunitsyn,  "Diffraction  Tomography  in  the  Fresnel  Zone", 

Moscow  University  Physics  Bulletin,  Vol.41,  No  2,  pp. 43-48,  1986. 

5.  V.E. Kunitsyn,  N.G. Preobrazhensky,  E.D. Tereshchenko, 

"Reconstruction  of  ionospheric  irregularities  structure  based 
on  the  radioprobing  data",  Dokl.Akad.Nauk  USSR,  Vol.306, 
pp. 575-579,  1989. 

6.  V.E. Kunitsyn,  E.D. Tereshchenko,  The  reconstruction  of  the 

ionosphere  Irregularities  structure ,  Preprint  Polar  Geophys. 

Inst . ,  J»  90-01  -69 ,  pp .  1  -60 ,  1 990 . 

7.  V.E. Kunitsyn,  "Diffraction  Tomography  based  on  small-angle 

scattering  data".  Transactions  SPIE,  Vol.1843,  pp. 172-182,  1992. 

8.  V.E. Kunitsyn,  E.D. Tereshchenko,  "Radio  Tomography  of  the 

Ionosphere",  Antennas  Rropagation  Magazine,  Vol.34,  pp. 22-32,  1992 

9.  V.L. Ginzburg,  Propagation  of  electromagnetic  waves  in 

plasma.  New  York,  Gordon  and  Breach  Science  Publishers,  Inc.,  1961 

10.  L.D. Landau,  E.M.Lifshits,  Statistical  physics.  New  York, 

Pergamon  Press,  1980. 

11.  L.D.Faddeev,  "Three-dimensional  inverse  scattering  problem", 

J.Sov.Math. ,  Vol  5,  pp.  334-360,  1976. 

12.  R.G. Newton,  "Inverse  scattering  II.  Three  dimentions"  J.Math. 

Phys.,  Vol. 21,  pp. 1698-1715,  1980. 

13.  K.Chadan,  P.C. Sabatier,  Inverse  problems  in  quantum  theory. 

New  York,  Springer-Verlag,  1977. 

14.  A. J.Devaney,  "Geophys.  diffraction  tomography",  IEEE  Trans. 

Geosci.  and  Remote  Sencing,  Vol. 22,  pp.3-13,  1984. 

15.  R.K.Mueller,  M.Kaveh,  G.Wade,  "Reconstructive  tomography  and 

applications  to  ultrasonics",  Proc.  IEEE,  Vol. 67,  pp. 567-587,  1979 

16.  D.Lesselier,  D.Vuillet-Laurent,  F.Jouvie  at  al.,  "Iterative 

solutions  of  some  direct  and  inverse  problems  in  electromagnetics 
and  acoustics".  Electromagnetics,  Vol. 5,  pp. 147-189,  1985. 


45" 


17.  K.  J.Langeberg,  "Applied,  inverse  problems  for  acoustic, 

electromagnetic  and  elastic  wave  scattering",  in  F. Sabatier  (ed. ), 
Basic  methods  of  tomography  and  inverse  problems,  Bristol, 

Adam  Hilger,  pp. 127-467,  1987. 

18.  W.Tabbara,  B.Duchene,  C.Pichot.  et  al,  "Diffraction  tomography: 

contribution  to  the  analysis  of  some  applications  in  microwaves 
and  ultrasonics".  Inverse  Problems,  Vol.4,  pp. 305-331 ,1988. 

19.  W.K. Pratt,  Digital  Image  Processing,  New  York,  A  Wiley 

Interscience  Publ.,  1978. 

20.  D.E.Dadgeon,  R.M.Mersereau,  Multidimensional  Digital  Signal 

Processing,  Prentice-Hall,  Inc.,  1984. 

21.  D.Slepian,  "About  bandwidth",  Proc.IEEE,  Vol.64,  No. 3,  pp. 292-300. 

22.  H.J. Landau,  H.O. Poliak,  "Prolate  spheroidal  wave  functions, 

Fourier  analysis  and  uncertainty".  Bell  Syst.Tech. J. ,  Vol.40, 
pp.  65-84;  Vol.41.  pp. 1295-1 336,  1961. 

23.  E.S. Andreeva,  V.E. Kunitsyn,  "The  reconstruction  of  the  projections 

of  scattered  irregularities  by  inexact  wavefield  data".  In: 

Waves  propagation  and  diffraction  in  inhomogeneous  media, 

Moscow,  Moscow  Inst.  Physics  and  Technology,  pp.  9-14,  1989. 

24.  V.E. Kunitsyn,  E.D. Tereshchenko,  "Determination  of  the 

turbulent  spectrum  in  the  Ionosphere  by  a  tomographic  method", 
Joum. Atm. Terr. Phy s. ,  Vol.54,  pp. 1275-1282,  1992. 

25.  A.V.Galinov,  V.E. Kunitsyn,  E.D. Tereshchenko,  "A  tomographic 

approach  to  the  investigation  of  a  randomly  nonuniform 
ionosphere".  Geomagnetism  and  Aeronomy,  Vol.31,  No. 3, 
pp. 351 -356,  1991 

26.  V.I.Tatarski,  Wave  propagation  in  a  turbulent  medium. , 

New  York,  McGraw-Hill,  1961. 

27.  A.Ishimaru,  Wave  propagation  and  scattering  in  random  media, 

Orlando  (Florida),  Academic  Press ;  1978. 

28.  I. M. Sobol,  Monte-Carlo  method,  Moscow,  Nauka,  1978  (in  Russian). 

29.  S.Chakrabartl,  S.K.Mitra,  "Calculations  of  two  dimensional 

digital  filters  by  spectral  convolution",  Proc.IEEE, 

Vol.65,  No. 6,  1977. 

30.  R. Gordon,  R. Bender,  A. Lent,  Joum. Theor. Biol. ,  Vol.29,  1970. 


A  (KM) 


Pc 

PL2 

0.154 

0.183 

0.345 

0.436 

0.439 

0.675 

0.672 

0.882 

0.093 

0.101 

0.241 

0.253 

0.468 

0.475 

0.654 

0.717 

5 

Calculation  of  realization 
of  power  spectrum. 

£(p.V 


r 


Fig.  2  3 


